MRS Longitudinal Modelling
Meaghan Perdue
12/14/2021
Organize data
Merge output data from multiple batches, select metabolite columns, add participant info and other data
acg_batch1 <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/ACG_LCM_CorrectedTable_12-14-2021-16-15.csv")
acg_batch2 <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/ACG_LCM_CorrectedTable_12-14-2021-16-21.csv")
acg_LCM_corr <- bind_rows(acg_batch1, acg_batch2)
anyDuplicated(acg_LCM_corr$File_ID) #check for duplicated values## [1] 106acg_LCM_corr<-unique(acg_LCM_corr) #remove duplicate rows
acg_LCM_corr["4", "Cho_GPC_PCh_conc_molar"]<-NA #Cho conc invalid for PS10008_PS14_046_2386 due to spike in spectrum
acg_LCM_corr["4", "Cho_GPC_PCh_conc_molal"]<-NA #Cho conc invalid for PS10008_PS14_046_2386 due to spike in spectrum
write.csv(acg_LCM_corr, "/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/ACG_LCM_CorrectedTable_12-14-2021_all.csv", row.names = FALSE)
acg_molal <- acg_LCM_corr %>%
dplyr::select(File_ID, fGM, fWM, fCSF, NAA_conc_molal, Cre_PCr_conc_molal, Cho_GPC_PCh_conc_molal, Ins_conc_molal, Glu_conc_molal, Glu_Gln_conc_molal)
acg_demo <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/acg_analysis_subs.csv")
acg_dat <- left_join(acg_demo, acg_molal, by = "File_ID")
write.csv(acg_dat, "/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/acg_data_molal_n113.csv", row.names = FALSE)
lag_batch1 <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/LAG_LCM_CorrectedTable_12-13-2021-15-18.csv")
lag_batch2 <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/LAG_LCM_CorrectedTable_12-14-2021-11-41.csv")
lag_batch3 <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/LAG_LCM_CorrectedTable_12-14-2021-12-18.csv")
lag_LCM_corr <- bind_rows(lag_batch1, lag_batch2, lag_batch3)
write.csv(lag_LCM_corr, "/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/LAG_LCM_CorrectedTable_12-14-2021_all.csv", row.names = FALSE)
lag_molal <- lag_LCM_corr %>%
dplyr::select(File_ID, fGM, fWM, fCSF, NAA_conc_molal, Cre_PCr_conc_molal, Cho_GPC_PCh_conc_molal, Ins_conc_molal, Glu_conc_molal, Glu_Gln_conc_molal)
lag_demo <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/lag_analysis_subs.csv")
lag_dat <- left_join(lag_demo, lag_molal, by = "File_ID")
lag_dat$subj_id<-as.character(lag_dat$subj_id) #make subj_id character type
write.csv(lag_dat, "/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/lag_data_molal_n321.csv", row.names = FALSE)##Check data distribution & outliers Create plots to check data distribution and identify outliers
ACG
plot_age <- ggplot(acg_dat, aes(mri_age_y))
plot_age + geom_density()
#age skewed
plot_naa <- ggplot(acg_dat, aes(NAA_conc_molal))
plot_naa + geom_density()
plot_naa + geom_boxplot()
rosnerTest(acg_dat$NAA_conc_molal, k = 4)## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4
## 2.759147 2.825357 2.694361 2.721635
##
## $sample.size
## [1] 113
##
## $parameters
## k
## 4
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4
## 3.425263 3.422302 3.419309 3.416284
##
## $n.outliers
## [1] 0
##
## $alternative
## [1] "Up to 4 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 13.36240 15.05227 14.57780 13.83339 14.71206 16.20111 15.44740 16.98446
## [9] 15.67072 15.95324 14.87793 12.74216 16.29728 14.71133 14.66048 13.28760
## [17] 14.59532 13.55773 16.85293 15.33699 12.45017 15.04043 15.27049 11.24177
## [25] 14.98401 14.18774 15.83415 13.91788 17.10858 12.09629 14.70242 14.40906
## [33] 16.46456 15.58993 15.03289 14.84971 11.30097 13.87568 15.24533 15.99486
## [41] 14.53574 12.16795 14.20751 14.18561 16.47436 14.52629 15.00648 16.20561
## [49] 13.22605 13.57437 14.05939 13.81423 17.33147 16.31303 16.29984 16.85100
## [57] 15.24641 14.93994 15.14031 15.13625 18.29604 13.89531 13.42919 15.16179
## [65] 13.87656 15.48760 12.64660 16.01636 13.64857 16.43270 15.59565 11.64927
## [73] 15.46789 14.18371 12.72496 16.73666 16.33378 15.15593 15.86369 15.28939
## [81] 15.31529 14.67805 13.38081 15.01442 13.01089 17.11696 14.54570 15.04389
## [89] 14.64002 15.57552 15.21758 15.80866 15.16132 14.92979 14.18939 16.96225
## [97] 13.86479 15.05179 15.22074 15.65731 15.51013 14.39392 15.03669 16.39137
## [105] 14.47850 14.65614 15.94476 16.78988 16.07591 13.54665 13.48426 14.73463
## [113] 15.68117
##
## $data.name
## [1] "acg_dat$NAA_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 14.88986 1.322180 11.24177 24 2.759147 3.425263 FALSE
## 2 1 14.92243 1.281772 11.30097 37 2.825357 3.422302 FALSE
## 3 2 14.95506 1.239990 18.29604 61 2.694361 3.419309 FALSE
## 4 3 14.92469 1.203472 11.64927 72 2.721635 3.416284 FALSE
##
## attr(,"class")
## [1] "gofOutlier"plot_cre <- ggplot(acg_dat, aes(Cre_PCr_conc_molal))
plot_cre + geom_density()
plot_cre + geom_boxplot()
rosnerTest(acg_dat$Cre_PCr_conc_molal, k = 5)## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4 R.5
## 3.394290 2.833607 2.707716 2.578893 2.649738
##
## $sample.size
## [1] 113
##
## $parameters
## k
## 5
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4 lambda.5
## 3.425263 3.422302 3.419309 3.416284 3.413225
##
## $n.outliers
## [1] 0
##
## $alternative
## [1] "Up to 5 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 10.869437 10.445542 14.138594 11.979460 11.842544 12.567296 12.183610
## [8] 13.315617 13.532807 13.910100 12.102060 9.609688 12.566391 11.689508
## [15] 12.650725 10.976876 11.888348 11.823341 12.975541 12.087174 10.745566
## [22] 12.005758 13.958950 11.012768 13.437536 12.359373 12.418194 10.775703
## [29] 13.584172 11.251205 11.993304 12.190118 14.756858 13.780875 12.273252
## [36] 11.243397 8.406874 11.279114 11.602557 13.761404 12.044563 10.276615
## [43] 10.918920 10.896575 11.867174 12.189638 12.580644 12.364280 11.517322
## [50] 10.442494 11.903800 12.566843 12.411654 12.733123 12.852863 14.066514
## [57] 10.742053 11.977272 12.534036 12.467809 13.778825 11.508522 11.740032
## [64] 11.569008 10.891687 12.474281 11.395342 12.716648 12.033736 13.646045
## [71] 12.873670 9.632426 11.901341 12.610355 9.187364 13.087309 12.160354
## [78] 11.717870 12.392428 11.411606 12.417766 11.214681 11.894967 11.689931
## [85] 11.069854 13.648131 12.837410 11.764882 12.329981 11.666319 12.825568
## [92] 13.032058 12.630243 12.355787 11.038718 9.731136 12.044142 11.736949
## [99] 11.772225 11.427900 12.588592 11.074387 11.520514 12.986814 12.563423
## [106] 12.334515 13.095183 12.401466 11.143206 11.414014 11.499070 12.331493
## [113] 11.900989
##
## $data.name
## [1] "acg_dat$Cre_PCr_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 12.03590 1.0691568 8.406874 37 3.394290 3.425263 FALSE
## 2 1 12.06830 1.0167044 9.187364 75 2.833607 3.422302 FALSE
## 3 2 12.09426 0.9833373 14.756858 33 2.707716 3.419309 FALSE
## 4 3 12.07005 0.9540396 9.609688 12 2.578893 3.416284 FALSE
## 5 4 12.09263 0.9284692 9.632426 72 2.649738 3.413225 FALSE
##
## attr(,"class")
## [1] "gofOutlier"plot_cho <- ggplot(acg_dat, aes(Cho_GPC_PCh_conc_molal))
plot_cho + geom_density()## Warning: Removed 1 rows containing non-finite values (stat_density).
plot_cho + geom_boxplot()## Warning: Removed 1 rows containing non-finite values (stat_boxplot).
rosnerTest(acg_dat$Cho_GPC_PCh_conc_molal, k = 4)## Warning in rosnerTest(acg_dat$Cho_GPC_PCh_conc_molal, k = 4): 1 observations
## with NA/NaN/Inf in 'x' removed.## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4
## 2.692059 2.486294 2.542571 2.587426
##
## $sample.size
## [1] 112
##
## $parameters
## k
## 4
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4
## 3.422302 3.419309 3.416284 3.413225
##
## $n.outliers
## [1] 0
##
## $alternative
## [1] "Up to 4 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 3.054728 2.338141 2.483788 2.615586 2.904744 3.029351 2.945605 2.946224
## [9] 2.469310 2.647417 2.329234 2.525650 2.964839 2.377024 2.316960 2.388015
## [17] 2.587590 2.964150 2.987201 1.884924 2.294818 2.160623 1.723047 2.273437
## [25] 2.486538 2.364455 2.110030 3.042283 2.145223 2.626855 2.657139 3.307897
## [33] 2.993478 2.528009 2.694001 2.751861 2.169544 2.576529 2.777803 2.298135
## [41] 2.093485 2.245856 1.887058 3.332061 2.683945 2.960541 2.306091 2.400287
## [49] 2.262987 2.373107 2.418340 2.428692 2.651362 2.724707 2.627568 3.213338
## [57] 2.541659 2.436647 2.503022 3.168533 2.273151 2.230455 2.575725 2.463493
## [65] 2.964749 2.427496 2.382196 2.398624 2.592440 2.075228 1.967159 2.210616
## [73] 2.178031 1.572637 2.265866 2.193475 2.435427 2.749229 2.738993 2.096637
## [81] 2.376787 2.150157 2.031053 1.736687 2.549619 2.032952 2.153716 2.354125
## [89] 2.585598 2.624180 2.597158 2.503254 2.679817 2.109408 2.755517 2.618999
## [97] 2.254582 2.379855 2.105926 2.389925 2.621861 2.349337 2.962557 2.620008
## [105] 2.643999 2.559238 3.322683 2.403398 2.301871 2.505683 2.762249 2.704160
##
## $data.name
## [1] "acg_dat$Cho_GPC_PCh_conc_molal"
##
## $bad.obs
## [1] 1
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 2.495888 0.3429536 1.572637 75 2.692059 3.422302 FALSE
## 2 1 2.504206 0.3329678 3.332061 45 2.486294 3.419309 FALSE
## 3 2 2.496680 0.3248691 3.322683 108 2.542571 3.416284 FALSE
## 4 3 2.489102 0.3164519 3.307897 33 2.587426 3.413225 FALSE
##
## attr(,"class")
## [1] "gofOutlier"plot_ins <- ggplot(acg_dat, aes(Ins_conc_molal))
plot_ins + geom_density()
shapiro.test(acg_dat$Ins_conc_molal) #not normal##
## Shapiro-Wilk normality test
##
## data: acg_dat$Ins_conc_molal
## W = 0.96211, p-value = 0.002754plot_ins + geom_boxplot()
rosnerTest(acg_dat$Ins_conc_molal, k = 3) #outlier value 1.471244, observation 84 (PS18_013) - spectrum showed issue, exclude from analysis ## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3
## 3.753416 3.044570 3.051178
##
## $sample.size
## [1] 113
##
## $parameters
## k
## 3
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3
## 3.425263 3.422302 3.419309
##
## $n.outliers
## [1] 1
##
## $alternative
## [1] "Up to 3 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 7.633881 8.178435 7.218240 7.696923 8.696692 8.600324 9.402207
## [8] 8.875539 6.779586 5.935284 10.304907 5.028532 7.058793 9.067489
## [15] 6.700483 7.361395 8.124760 8.682546 6.553380 6.320046 7.315123
## [22] 7.446418 7.281237 7.102747 9.505161 7.839894 7.887706 7.616099
## [29] 7.601724 6.347168 7.912922 5.656884 9.656654 9.672370 8.341248
## [36] 8.816154 5.196891 6.304340 6.789199 8.018694 7.594973 7.416233
## [43] 7.406756 3.198494 2.985922 9.898020 8.079535 8.703354 7.157223
## [50] 6.213808 8.399767 9.174211 9.656691 9.335033 10.080614 4.813291
## [57] 9.728369 7.521052 7.441993 8.051283 10.887902 8.695335 5.872553
## [64] 6.418561 5.242212 8.749670 7.765345 8.185832 7.950928 8.394652
## [71] 9.583041 6.202420 7.829477 6.209769 4.285017 8.677590 7.775045
## [78] 9.334317 8.684385 7.236647 7.742767 6.219638 9.105159 1.471244
## [85] 5.218826 9.407651 7.634559 6.064850 8.785683 3.770173 7.618745
## [92] 8.763466 7.123140 6.399259 6.804541 5.855502 9.364210 6.203701
## [99] 10.565682 9.856594 7.453487 7.728481 7.846516 9.947994 8.640216
## [106] 7.758280 8.476271 8.672477 6.437713 5.739577 6.673856 7.596796
## [113] 7.229168
##
## $data.name
## [1] "acg_dat$Ins_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 7.588864 1.629881 1.471244 84 3.753416 3.425263 TRUE
## 2 1 7.643485 1.529793 2.985922 45 3.044570 3.422302 FALSE
## 3 2 7.685445 1.470564 3.198494 44 3.051178 3.419309 FALSE
##
## attr(,"class")
## [1] "gofOutlier"acg_dat["84", "Ins_conc_molal"] <- NA #replace outlier value with NA to exclude from analysis
plot_glx <- ggplot(acg_dat, aes(Glu_Gln_conc_molal))
plot_glx + geom_density()
plot_glx + geom_boxplot()
rosnerTest(acg_dat$Glu_Gln_conc_molal, k = 1)## $distribution
## [1] "Normal"
##
## $statistic
## R.1
## 2.840029
##
## $sample.size
## [1] 113
##
## $parameters
## k
## 1
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1
## 3.425263
##
## $n.outliers
## [1] 0
##
## $alternative
## [1] "Up to 1 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 22.16208 24.54895 30.79902 25.15445 23.12144 23.45996 27.95219 28.85494
## [9] 28.47293 30.83024 24.02170 20.83895 26.81225 25.81936 25.83210 28.15351
## [17] 28.69984 24.45281 25.87441 26.10568 24.90431 27.69072 29.16664 25.73952
## [25] 31.00567 25.29233 29.35298 28.81394 30.45098 26.42377 25.25458 30.67838
## [33] 29.65833 34.25604 26.24353 28.06087 21.93668 25.00839 27.05607 29.92424
## [41] 27.12366 23.65348 26.82576 24.86136 26.92752 23.98651 26.21693 28.06806
## [49] 26.26637 27.44238 23.98526 25.40089 22.59948 23.46421 28.23542 33.70878
## [57] 26.99835 25.59893 26.95525 25.77541 32.34599 24.66776 25.46001 25.71687
## [65] 29.73266 25.54233 27.95401 29.18637 26.09421 27.74434 27.70571 21.88957
## [73] 26.72328 30.44184 21.32831 31.44932 28.73199 25.08766 29.57015 24.10930
## [81] 27.84161 24.85586 29.21009 20.93551 23.34771 27.38610 26.35670 26.10330
## [89] 28.90132 26.67329 24.94565 29.44010 29.32746 26.23233 24.32566 24.39736
## [97] 26.88103 30.91320 28.56651 27.47953 25.89977 24.22166 25.89081 28.12079
## [105] 29.77257 28.04112 29.21096 29.61668 27.08694 28.10033 24.05989 26.66364
## [113] 26.45077
##
## $data.name
## [1] "acg_dat$Glu_Gln_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 26.82935 2.615006 34.25604 34 2.840029 3.425263 FALSE
##
## attr(,"class")
## [1] "gofOutlier"write.csv(acg_dat, "/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/acg_data_molal_outliers_removed_n113.csv", row.names = FALSE)LAG
plot_age <- ggplot(lag_dat, aes(mri_age_y))
plot_age + geom_density()
#age skewed
plot_naa <- ggplot(lag_dat, aes(NAA_conc_molal))
plot_naa + geom_density()
plot_naa + geom_boxplot()
rosnerTest(lag_dat$NAA_conc_molal, k = 10) #outliers val. 19.207, obs. 296 & Val. 10.148, obs. 306 - some issues in spectra, exclude from analysis for NAA## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4 R.5 R.6 R.7 R.8
## 4.467031 4.044357 3.599759 3.613937 3.039885 2.720506 2.647000 2.651977
## R.9 R.10
## 2.660497 2.662447
##
## $sample.size
## [1] 321
##
## $parameters
## k
## 10
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4 lambda.5 lambda.6 lambda.7 lambda.8
## 3.742579 3.741706 3.740829 3.739949 3.739067 3.738181 3.737291 3.736399
## lambda.9 lambda.10
## 3.735503 3.734604
##
## $n.outliers
## [1] 2
##
## $alternative
## [1] "Up to 10 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 14.16473 13.06076 14.35969 13.27563 16.92063 13.92045 14.90883 14.96705
## [9] 13.87538 13.79928 15.43119 13.74363 13.42774 14.32941 15.29816 15.36852
## [17] 14.96343 14.01462 15.11779 16.29892 15.28430 15.17137 15.34158 15.28575
## [25] 13.61378 16.76759 13.63702 14.41263 13.69061 15.64950 14.60942 15.25717
## [33] 15.68290 14.24180 15.96102 13.19084 14.22355 14.45688 13.89590 16.68967
## [41] 15.14115 13.88211 14.56337 14.77687 14.00542 14.12098 13.36743 13.70315
## [49] 14.59050 15.55276 13.33791 14.13365 12.08635 14.08027 13.90090 13.63249
## [57] 15.72428 11.86871 14.17563 11.96658 12.37485 14.57472 13.54554 14.35760
## [65] 13.86325 16.89273 14.68644 16.67785 15.93255 16.85334 14.93659 14.21829
## [73] 13.99673 14.23697 13.86957 13.86406 14.36716 14.41372 12.45632 14.76312
## [81] 14.75832 14.92673 13.88514 14.24089 15.42366 14.86462 14.32717 14.96667
## [89] 13.72426 14.16948 15.16894 13.93828 12.92563 13.94977 13.64160 14.55102
## [97] 14.51660 12.52392 14.21054 13.93670 16.19979 14.45536 14.10612 13.35969
## [105] 16.58862 14.41649 14.80811 14.26277 12.37256 15.44893 13.43728 13.57859
## [113] 13.72929 14.86513 11.95501 14.39891 14.11323 14.14168 14.73091 14.40429
## [121] 13.19741 13.71443 14.26528 14.92876 14.67063 14.35362 15.12298 13.95432
## [129] 13.68703 13.44733 14.01869 14.75568 14.74986 11.38458 12.08162 11.74260
## [137] 13.93366 14.16181 14.38636 14.40623 14.15278 14.17590 15.46643 14.09616
## [145] 13.62893 14.14878 13.54441 13.09883 12.79401 14.64345 14.02714 14.59131
## [153] 14.51512 15.41640 13.41648 13.54853 13.78080 13.79989 13.51420 14.77851
## [161] 15.60295 13.13390 13.99763 14.27699 15.00866 15.01634 12.09852 14.98452
## [169] 13.60870 16.34165 11.91384 12.33377 14.09992 15.06924 13.87378 13.47019
## [177] 14.02295 18.06089 14.25601 14.37746 13.24330 15.14502 13.95585 13.43714
## [185] 14.30145 16.30323 13.44700 14.70962 13.78765 13.43206 14.86922 14.73623
## [193] 14.54106 13.30455 14.61897 13.91250 15.44818 13.00756 14.21959 16.58017
## [201] 14.30707 15.14338 12.83493 13.96649 15.41109 13.77502 12.24289 15.19071
## [209] 13.32940 13.11120 15.00722 15.15742 13.96143 14.70099 14.27321 15.14983
## [217] 14.79365 15.19557 14.84745 13.19832 15.12751 13.95088 13.28456 14.90679
## [225] 12.94820 13.28282 14.39416 13.76430 14.91352 15.55618 15.20787 14.89958
## [233] 14.31622 14.53646 15.29391 13.51980 15.05685 14.73951 15.01094 15.61748
## [241] 14.09588 14.63147 15.59956 15.33848 14.07436 14.30913 13.94217 16.47833
## [249] 15.29907 13.85813 15.73137 15.08295 16.15169 14.28221 15.30824 14.47596
## [257] 14.47137 14.28477 15.33016 14.75400 14.54703 14.66072 13.49926 14.47973
## [265] 13.83620 14.83537 15.05792 14.18100 13.68113 13.28288 14.53224 13.92267
## [273] 14.80805 13.25470 14.60733 13.77756 15.36693 16.13917 14.36568 14.54767
## [281] 14.51644 16.31083 15.95003 13.41067 15.50889 13.79439 14.53462 14.83200
## [289] 14.42306 13.93431 13.79584 14.91428 12.66542 14.98849 15.38162 19.20747
## [297] 13.74790 17.99338 14.95991 13.93438 15.42154 13.27326 13.14415 14.47594
## [305] 13.99067 10.14812 13.73230 14.91913 15.85376 13.23125 13.40927 14.90803
## [313] 14.30314 14.54287 14.33013 16.55184 14.94880 14.30460 15.31501 14.46473
## [321] 16.31675
##
## $data.name
## [1] "lag_dat$NAA_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 14.39149 1.0781156 19.20747 296 4.467031 3.742579 TRUE
## 2 1 14.37644 1.0454867 10.14812 306 4.044357 3.741706 TRUE
## 3 2 14.38970 1.0198436 18.06089 178 3.599759 3.740829 FALSE
## 4 3 14.37815 1.0003560 17.99338 298 3.613937 3.739949 FALSE
## 5 4 14.36675 0.9810139 11.38458 134 3.039885 3.739067 FALSE
## 6 5 14.37619 0.9680502 11.74260 136 2.720506 3.738181 FALSE
## 7 6 14.38455 0.9580956 16.92063 5 2.647000 3.737291 FALSE
## 8 7 14.37647 0.9488234 16.89273 66 2.651977 3.736399 FALSE
## 9 8 14.36843 0.9395707 11.86871 58 2.660497 3.735503 FALSE
## 10 9 14.37644 0.9303092 16.85334 70 2.662447 3.734604 FALSE
##
## attr(,"class")
## [1] "gofOutlier"lag_dat["296", "NAA_conc_molal"]<- NA
lag_dat["306", "NAA_conc_molal"]<- NA
plot_cre <- ggplot(lag_dat, aes(Cre_PCr_conc_molal))
plot_cre + geom_density()
plot_cre + geom_boxplot()
rosnerTest(lag_dat$Cre_PCr_conc_molal, k = 10) #outliers val. 14.829, obs. 241 (PS17_047) - structure in resid. & val. 14.09, obs. 178 (PS17_045, spectrum looked fine), excluding only PS17_047 for now## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4 R.5 R.6 R.7 R.8
## 4.835902 4.220654 3.332561 3.278883 3.286432 3.307998 2.911341 2.895521
## R.9 R.10
## 2.767415 2.775535
##
## $sample.size
## [1] 321
##
## $parameters
## k
## 10
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4 lambda.5 lambda.6 lambda.7 lambda.8
## 3.742579 3.741706 3.740829 3.739949 3.739067 3.738181 3.737291 3.736399
## lambda.9 lambda.10
## 3.735503 3.734604
##
## $n.outliers
## [1] 2
##
## $alternative
## [1] "Up to 10 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 9.979902 8.752796 9.874363 8.961297 11.510180 10.121568 10.418593
## [8] 11.028922 10.171413 9.929791 10.871830 11.064228 10.111309 9.319128
## [15] 10.662947 10.890644 10.197139 9.919287 10.426231 11.325982 10.557574
## [22] 11.088082 11.039646 11.172931 8.936929 12.109674 10.000368 10.420059
## [29] 9.572731 10.511121 9.984628 10.695378 11.039338 9.822119 11.097416
## [36] 9.239262 9.742404 10.192047 10.845538 12.117826 10.389433 9.535353
## [43] 10.010295 10.823639 9.876324 10.019161 8.608303 9.442680 11.444015
## [50] 11.419644 9.408142 10.139106 8.483950 9.792171 10.032519 8.968350
## [57] 11.562417 9.156420 10.285179 7.449145 10.159106 10.623228 10.245682
## [64] 11.514351 10.182653 11.786704 11.560836 12.034078 11.166365 11.463697
## [71] 10.397742 10.176308 11.796639 10.733740 10.174742 10.465407 9.824296
## [78] 11.420357 7.832529 10.284534 12.513826 11.576531 10.541650 10.808215
## [85] 10.483803 10.041962 10.382300 9.976776 10.765842 10.520082 11.422196
## [92] 10.843472 9.968723 10.781936 9.841352 9.930816 11.494602 9.908886
## [99] 10.765624 10.122441 11.914466 10.403873 10.056177 9.179081 11.252406
## [106] 10.011672 10.666794 9.981113 9.225683 10.367639 10.879950 11.335061
## [113] 9.425501 10.725975 10.241112 9.317008 9.889557 9.560765 10.593153
## [120] 11.263291 8.935236 10.137273 10.048855 11.211108 10.481657 10.028783
## [127] 11.460008 9.844337 10.207849 8.876822 9.557376 9.418484 9.631575
## [134] 11.004550 7.414051 8.984416 9.798072 10.236338 10.150536 10.108983
## [141] 10.653763 10.137920 10.541536 9.525275 9.886573 9.896683 9.185657
## [148] 9.642602 7.281217 9.968306 10.715037 10.739283 11.039381 11.250002
## [155] 9.759444 9.643949 9.675426 9.621214 9.141950 10.356422 10.788507
## [162] 9.362209 9.822128 9.645301 10.157630 10.381588 8.913720 11.109833
## [169] 9.890247 12.488886 9.115647 9.814255 10.627685 11.734993 10.295594
## [176] 10.173568 9.988731 14.090723 10.829617 10.247883 10.533434 10.999808
## [183] 10.001382 9.229915 10.354994 11.296388 9.994843 10.904378 9.345731
## [190] 9.618365 9.759703 10.530182 10.844392 8.895843 10.110129 11.423088
## [197] 11.096328 9.958072 9.931007 10.921273 9.154313 10.482202 10.423407
## [204] 10.629254 12.212044 9.913344 8.305599 10.446897 9.841379 9.379922
## [211] 9.844787 10.968766 10.801476 9.891839 10.632212 9.671280 10.872873
## [218] 10.666877 10.645958 9.493899 13.103099 9.512697 10.276565 9.846960
## [225] 9.773848 9.929345 11.034653 10.095542 10.422810 10.640722 10.407355
## [232] 11.450386 10.284795 10.179735 11.263605 10.183927 10.257883 10.937157
## [239] 10.529431 11.012918 14.829071 9.133491 11.531961 10.173887 9.610385
## [246] 9.349643 9.324653 10.926724 11.044735 8.859819 11.228731 10.816240
## [253] 11.271292 9.861827 10.852278 10.764773 9.380198 9.481786 10.781017
## [260] 8.985104 9.888393 9.434081 8.870117 9.122209 9.542512 9.889146
## [267] 10.769356 9.616033 8.966342 9.255494 10.764599 9.897261 9.427714
## [274] 8.163350 10.017130 9.265501 10.944724 12.004282 9.764853 9.856451
## [281] 10.020483 12.009948 10.044769 9.667192 11.197137 9.964094 10.035066
## [288] 10.812708 9.728669 9.918824 8.843802 10.202045 9.183452 9.822267
## [295] 10.790509 12.454846 8.641534 9.507122 10.361671 9.646740 11.011202
## [302] 9.075410 10.166484 9.737467 9.615290 7.882414 10.377305 10.224366
## [309] 12.449600 10.083783 9.260131 11.009901 9.914846 9.659798 10.610323
## [316] 11.561746 11.077058 9.291690 10.160005 9.934510 11.934339
##
## $data.name
## [1] "lag_dat$Cre_PCr_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 10.26105 0.9446048 14.829071 241 4.835902 3.742579 TRUE
## 2 1 10.24678 0.9107457 14.090723 178 4.220654 3.741706 TRUE
## 3 2 10.23473 0.8862592 7.281217 149 3.332561 3.740829 FALSE
## 4 3 10.24402 0.8719682 13.103099 221 3.278883 3.739949 FALSE
## 5 4 10.23500 0.8583617 7.414051 135 3.286432 3.739067 FALSE
## 6 5 10.24393 0.8448557 7.449145 60 3.307998 3.738181 FALSE
## 7 6 10.25280 0.8313244 7.832529 79 2.911341 3.737291 FALSE
## 8 7 10.26051 0.8213002 7.882414 306 2.895521 3.736399 FALSE
## 9 8 10.26810 0.8114875 12.513826 81 2.767415 3.735503 FALSE
## 10 9 10.26091 0.8027211 12.488886 170 2.775535 3.734604 FALSE
##
## attr(,"class")
## [1] "gofOutlier"lag_dat["241", "Cre_PCr_conc_molal"]<-NA #replace outlier value with NA to exclude from analysis
plot_cho <- ggplot(lag_dat, aes(Cho_GPC_PCh_conc_molal))
plot_cho + geom_density()
plot_cho + geom_boxplot()
rosnerTest(lag_dat$Cho_GPC_PCh_conc_molal, k = 8)## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4 R.5 R.6 R.7 R.8
## 3.694967 3.121637 3.061660 3.062209 2.979166 3.001536 2.871646 2.733239
##
## $sample.size
## [1] 321
##
## $parameters
## k
## 8
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4 lambda.5 lambda.6 lambda.7 lambda.8
## 3.742579 3.741706 3.740829 3.739949 3.739067 3.738181 3.737291 3.736399
##
## $n.outliers
## [1] 0
##
## $alternative
## [1] "Up to 8 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 2.068181 2.062467 2.386391 2.355035 2.498191 2.256110 2.077815 2.425292
## [9] 2.398770 1.895982 1.865260 2.184158 2.004899 2.128662 2.236830 2.378977
## [17] 2.355458 2.266506 2.404117 2.307009 2.052566 2.364941 2.267604 1.940308
## [25] 2.587651 2.716570 2.469831 2.547177 2.087747 2.170207 2.232363 2.241430
## [33] 2.441839 2.264953 2.333662 2.359801 2.270399 2.241648 1.920730 2.467142
## [41] 2.242001 2.227010 2.229581 2.053798 1.939191 2.364334 2.678152 2.362113
## [49] 1.676867 2.166236 2.306276 2.211050 2.178445 2.361600 2.344510 2.268087
## [57] 2.667157 2.059386 2.034573 2.045022 2.214431 2.995506 1.846508 2.183946
## [65] 2.072392 2.429630 2.330074 2.491408 2.265040 2.154392 2.411001 2.046897
## [73] 2.534412 2.027012 2.376839 2.222266 2.306644 2.292419 2.108096 2.172672
## [81] 1.839589 1.859195 1.443755 1.640698 2.420001 2.337525 2.238790 2.240345
## [89] 2.290795 2.399979 2.350376 2.151479 2.142910 2.110996 1.985190 1.928934
## [97] 2.304258 1.765061 2.137595 1.985591 2.781876 2.123475 1.756957 1.700318
## [105] 1.975316 1.974756 1.999963 2.031005 1.687589 1.708882 2.296279 2.182991
## [113] 2.276755 2.101586 1.837847 2.371711 2.065515 2.499424 2.363551 2.223515
## [121] 2.232809 2.045991 1.947255 2.339091 2.193248 2.528316 2.691920 2.777161
## [129] 2.451074 2.763753 2.377191 2.160428 2.215103 2.214312 2.533792 2.384939
## [137] 2.179386 2.552211 2.673232 2.476304 1.625876 1.926042 2.202001 2.204096
## [145] 2.344069 2.152159 1.895155 1.710908 1.401447 1.946252 2.047116 2.456831
## [153] 2.086612 2.201666 2.190181 2.294616 2.581968 2.225911 2.498999 2.000618
## [161] 2.163682 1.960245 2.204252 2.159811 2.097381 2.254865 1.844443 1.598729
## [169] 1.951488 3.208685 1.836756 1.950559 1.992121 2.281328 2.008449 2.599790
## [177] 2.317712 3.033330 2.698335 2.332814 2.583764 2.234143 2.191218 2.264957
## [185] 1.977277 2.067888 2.356329 2.431505 2.779803 2.150897 2.566469 2.101205
## [193] 2.159059 2.164490 2.152499 2.489624 2.791537 1.673365 2.193421 2.761471
## [201] 1.959300 2.022674 1.970196 2.366537 2.774857 2.129388 1.836375 2.360163
## [209] 2.476086 2.665161 2.258604 2.298965 2.012883 2.201029 2.158646 2.353116
## [217] 1.946114 1.912304 1.692442 2.307008 2.449170 2.132806 2.045527 2.155850
## [225] 2.166616 2.440956 1.783126 2.272129 2.595674 2.284913 2.155251 2.400964
## [233] 2.048543 2.045160 2.356748 2.222524 2.346513 2.378440 2.281038 1.906127
## [241] 1.725583 2.479767 2.262261 2.082256 2.047599 2.138599 2.252565 2.262486
## [249] 1.881557 2.322610 1.489680 2.225178 2.118903 2.390375 2.437536 2.400945
## [257] 1.832693 2.290250 2.007843 2.346441 2.089229 2.428740 2.391081 2.732220
## [265] 1.971298 2.235255 2.116196 1.941129 2.009743 1.999019 1.927678 1.994412
## [273] 2.372282 2.491822 1.886355 2.109659 2.028203 1.807961 2.439615 2.152157
## [281] 1.892138 2.163111 2.376896 1.914835 2.734225 1.944843 1.756721 2.135482
## [289] 2.101386 1.768588 1.974516 1.975407 2.167757 2.112266 2.510332 2.521732
## [297] 2.431779 2.959724 2.046249 1.986009 2.300110 2.368720 2.151772 2.070465
## [305] 2.457353 2.243696 2.105342 2.108232 2.180180 2.349512 1.911659 2.150353
## [313] 2.011103 2.243254 2.417149 1.800195 2.190437 2.875345 2.116343 2.319437
## [321] 2.422916
##
## $data.name
## [1] "lag_dat$Cho_GPC_PCh_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 2.208213 0.2707661 3.208685 170 3.694967 3.742579 FALSE
## 2 1 2.205087 0.2653235 3.033330 178 3.121637 3.741706 FALSE
## 3 2 2.202490 0.2616371 1.401447 149 3.061660 3.740829 FALSE
## 4 3 2.205009 0.2581460 2.995506 62 3.062209 3.739949 FALSE
## 5 4 2.202516 0.2546890 1.443755 83 2.979166 3.739067 FALSE
## 6 5 2.204917 0.2514736 2.959724 298 3.001536 3.738181 FALSE
## 7 6 2.202521 0.2482341 1.489680 251 2.871646 3.737291 FALSE
## 8 7 2.204791 0.2453333 2.875345 318 2.733239 3.736399 FALSE
##
## attr(,"class")
## [1] "gofOutlier"plot_ins <- ggplot(lag_dat, aes(Ins_conc_molal))
plot_ins + geom_density()
plot_ins + geom_boxplot()
rosnerTest(lag_dat$Ins_conc_molal, k = 8) #outliers: val. 14.077, obs. 241 (PS17_047) - structure in resid; val. 11.777, obs. 221 (PS16_014) - spectrum fine, keep for now; val. 1.542, obs. 197 (PS18_028) - noisy spectrum, spiky residuals, exclude; val. 2.04, obs. 62 (PS14_111) - spectrum fine, keep for now## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4 R.5 R.6 R.7 R.8
## 5.981441 4.463363 4.077869 3.761599 3.384729 3.303295 3.148162 3.105675
##
## $sample.size
## [1] 321
##
## $parameters
## k
## 8
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4 lambda.5 lambda.6 lambda.7 lambda.8
## 3.742579 3.741706 3.740829 3.739949 3.739067 3.738181 3.737291 3.736399
##
## $n.outliers
## [1] 4
##
## $alternative
## [1] "Up to 8 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 8.084411 7.069186 6.913031 6.915494 7.726417 8.046174 6.937152
## [8] 4.793633 5.952276 6.984233 6.400936 4.413732 7.176021 7.448535
## [15] 7.547442 6.541306 8.400250 5.957678 6.965664 7.688786 6.245160
## [22] 7.360921 7.925104 4.710055 6.035369 7.937845 7.166536 7.269577
## [29] 5.736145 7.751891 6.217281 7.121735 7.138640 6.104017 7.736566
## [36] 6.539681 4.958196 7.275083 6.485777 6.598573 7.153802 7.352086
## [43] 7.119035 6.146843 6.038556 6.638197 4.629802 6.187281 6.553798
## [50] 5.945866 7.474194 5.559932 5.743412 6.163849 8.347824 4.432658
## [57] 8.169688 5.045280 8.873237 6.213139 6.400722 2.040275 8.158702
## [64] 6.760986 6.468052 7.273556 6.742605 7.552488 7.843304 5.876394
## [71] 8.162408 4.963944 7.945325 5.440451 4.854219 4.834565 5.839066
## [78] 5.839637 5.890621 5.871859 6.963936 5.579543 5.190653 7.062975
## [85] 6.065434 7.046003 6.865074 6.210505 6.024185 6.121049 8.133551
## [92] 6.863480 6.537634 7.457231 5.008241 6.892333 5.936366 3.826364
## [99] 7.075620 5.567096 5.248858 5.491554 4.529573 7.011630 5.896171
## [106] 6.268198 5.292489 5.663161 6.753982 5.643647 6.788746 5.437125
## [113] 3.078293 7.511207 2.568354 8.286250 5.309958 6.053735 8.092627
## [120] 7.122038 6.509101 6.564424 5.780997 6.770400 5.727577 7.829891
## [127] 6.206182 7.258719 6.558572 4.787960 6.132452 6.029208 6.663446
## [134] 7.630692 4.608510 5.731585 5.435869 5.209567 2.734199 5.857696
## [141] 5.638088 6.135153 6.888234 5.236533 5.655180 4.948126 4.992038
## [148] 6.025793 2.972647 6.870427 7.728083 7.800580 7.334580 6.676194
## [155] 6.170230 5.196197 7.073851 6.127912 7.232895 4.960469 6.612268
## [162] 6.088183 5.798318 6.543746 6.997835 6.987522 4.250557 4.637685
## [169] 6.778603 7.502499 4.679992 5.270903 6.367636 6.301188 7.918660
## [176] 6.533953 7.497948 8.801777 6.797100 5.264895 5.014240 6.141243
## [183] 5.472028 3.885218 4.816043 6.248540 4.626087 6.007857 6.538200
## [190] 7.586651 7.411192 7.320786 7.682231 7.120925 6.686811 7.373729
## [197] 1.541994 6.314608 6.840831 7.095049 6.617076 6.056406 5.882616
## [204] 4.775951 5.731427 5.979244 6.418407 6.735628 5.243289 5.869558
## [211] 5.226573 6.275003 7.093997 7.541310 5.731832 4.649731 6.150065
## [218] 6.428658 7.426624 6.729085 11.777066 5.591387 3.080430 8.027587
## [225] 7.066068 6.661947 6.600256 6.955873 6.949031 6.824263 7.605754
## [232] 6.414457 6.094616 6.805794 5.865261 8.270934 4.935241 6.818562
## [239] 5.441384 7.674150 14.076805 4.939637 6.137708 5.300423 6.306576
## [246] 6.245176 7.663984 7.103391 8.335725 3.649002 9.145931 7.122485
## [253] 8.331023 7.199990 8.067478 7.557525 5.521439 6.036198 6.082461
## [260] 7.524328 5.926927 5.705055 5.399996 6.232379 5.707127 4.605042
## [267] 4.333050 5.607343 5.940738 5.740379 6.177392 6.685727 5.210419
## [274] 6.865329 6.749236 5.464544 6.496419 6.705573 7.274507 6.221882
## [281] 5.994352 6.059500 7.660503 6.660801 7.439851 6.777698 5.282780
## [288] 8.773125 6.173129 5.673978 5.639679 5.174471 8.615918 5.774767
## [295] 5.411721 7.832813 7.925781 4.996084 6.276000 6.161567 6.484734
## [302] 5.594141 6.793109 6.401047 9.077790 8.562826 5.674753 6.368745
## [309] 8.425749 6.658481 4.433254 5.482190 7.149805 6.079902 6.139307
## [316] 7.111624 5.971685 6.120964 6.571466 5.418973 6.379001
##
## $data.name
## [1] "lag_dat$Ins_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 6.383077 1.286267 14.076805 241 5.981441 3.742579 TRUE
## 2 1 6.359034 1.213890 11.777066 221 4.463363 3.741706 TRUE
## 3 2 6.342050 1.177099 1.541994 197 4.077869 3.740829 TRUE
## 4 3 6.357144 1.147615 2.040275 62 3.761599 3.739949 TRUE
## 5 4 6.370762 1.123401 2.568354 115 3.384729 3.739067 FALSE
## 6 5 6.382795 1.104533 2.734199 139 3.303295 3.738181 FALSE
## 7 6 6.394378 1.086898 2.972647 149 3.148162 3.737291 FALSE
## 8 7 6.405275 1.071259 3.078293 113 3.105675 3.736399 FALSE
##
## attr(,"class")
## [1] "gofOutlier"lag_dat["241", "Ins_conc_molal"] <- NA #replace outlier value with NA to exclude from analysis
lag_dat["197", "Ins_conc_molal"] <- NA
plot_glx <- ggplot(lag_dat, aes(Glu_Gln_conc_molal))
plot_glx + geom_density()
plot_glx + geom_boxplot()
rosnerTest(lag_dat$Glu_Gln_conc_molal, k = 7)#outliers: val. 10.09, obs. 241 (PS17_047) - structure in resid; val. 31.2899, obs. 296 (PS21_030) - big spike in resid. around 2.1, exclude.## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4 R.5 R.6 R.7
## 4.954047 4.247538 3.396313 3.439741 3.463175 3.267721 3.232088
##
## $sample.size
## [1] 321
##
## $parameters
## k
## 7
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4 lambda.5 lambda.6 lambda.7
## 3.742579 3.741706 3.740829 3.739949 3.739067 3.738181 3.737291
##
## $n.outliers
## [1] 2
##
## $alternative
## [1] "Up to 7 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 23.58311 20.72833 21.43465 20.79838 24.55542 19.38008 22.13533 21.98148
## [9] 18.58332 21.80094 23.39061 22.39400 22.25994 19.84839 22.62240 19.74196
## [17] 19.55718 20.18644 20.93614 24.01329 22.14911 23.96545 23.49753 23.68093
## [25] 20.83139 25.55075 18.83744 23.02848 20.10324 24.84118 23.11824 25.62522
## [33] 25.14532 20.60338 23.10183 18.68727 19.90887 22.63726 23.44172 23.97287
## [41] 20.79684 21.03860 21.41871 22.38334 19.07457 19.80466 21.54623 21.54133
## [49] 22.80842 25.04145 20.65555 22.02043 18.52301 21.98565 22.88652 20.99602
## [57] 23.64194 20.45706 24.44658 22.40116 19.22839 24.77315 20.40702 22.32472
## [65] 22.79321 23.27919 23.66026 23.59194 23.28915 23.15967 20.72270 22.14974
## [73] 23.87191 27.09158 23.12072 25.11318 22.08527 21.17324 14.25503 21.44083
## [81] 23.08444 21.47161 23.79592 26.39922 24.52352 23.18206 22.57951 22.53444
## [89] 24.24539 23.38651 21.39081 23.50524 21.84539 23.16617 20.66627 20.68771
## [97] 21.33913 23.82111 24.17549 20.43359 24.17777 21.80255 21.79162 21.34246
## [105] 26.43216 21.02534 21.61161 22.57479 21.58113 22.47585 21.89183 22.86451
## [113] 21.63056 23.83333 23.16833 22.80824 21.74170 19.79860 22.68433 25.43932
## [121] 20.46600 19.55287 21.71921 23.91575 19.94617 22.54020 25.12373 20.26601
## [129] 20.42848 18.40115 22.31740 21.34145 21.18639 24.46723 17.55526 19.17975
## [137] 22.14220 21.22080 20.62138 19.29726 22.33080 22.36491 23.64798 24.04583
## [145] 25.02751 22.02316 20.84404 22.82850 23.75262 19.93399 21.93649 21.46368
## [153] 23.51970 23.79151 22.29261 22.55887 19.88319 23.52534 18.26039 22.44651
## [161] 23.87196 21.33171 21.19221 21.21896 25.13735 21.72890 21.00404 24.38403
## [169] 22.77448 24.07955 20.92792 22.79597 21.53526 25.58587 20.23221 20.38536
## [177] 19.64011 29.12643 21.85437 22.81675 24.27197 20.09271 21.28381 21.20096
## [185] 20.61092 26.55674 20.24095 22.22341 21.17598 21.05662 20.70603 22.17497
## [193] 22.22738 17.77073 21.62541 23.02400 23.16615 17.38515 20.64300 21.97048
## [201] 23.41849 22.42124 20.67899 21.35918 23.92414 18.75710 19.72431 22.35249
## [209] 23.13971 19.66645 21.97753 21.86935 22.15373 18.04793 20.63749 19.81250
## [217] 22.91987 21.99328 21.86022 19.23943 22.18618 21.97414 22.10059 21.96999
## [225] 21.38803 21.08812 21.12026 20.43893 20.01495 22.48505 23.44854 23.27656
## [233] 22.72471 21.03119 22.44160 19.23947 24.24070 21.72010 21.19354 21.24039
## [241] 10.09497 19.32810 23.68178 22.11865 18.91616 19.59606 19.02554 22.92979
## [249] 26.99382 18.91880 22.36158 21.71252 22.51005 20.90044 23.69281 22.60270
## [257] 20.69658 19.00693 19.97580 15.08181 17.70438 17.08767 17.40514 16.78328
## [265] 20.40818 19.57305 21.57031 21.73287 18.10301 16.64465 21.45075 20.54401
## [273] 19.32920 14.36555 21.66594 17.20659 21.57909 24.57946 17.49328 18.35123
## [281] 21.22156 24.11680 22.51873 19.19853 23.08553 19.03798 23.38476 21.59376
## [289] 19.09657 21.28636 18.16122 19.70762 18.41625 18.73419 22.65583 31.28994
## [297] 18.89292 17.32671 18.83142 19.60466 24.66522 18.35557 22.92964 18.49429
## [305] 21.56271 18.98307 20.17174 24.36413 25.12361 21.76663 19.89145 21.01574
## [313] 20.06260 25.41609 26.15708 23.69646 23.29001 20.54746 23.65148 20.44966
## [321] 28.51150
##
## $data.name
## [1] "lag_dat$Glu_Gln_conc_molal"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 21.68857 2.340229 10.09497 241 4.954047 3.742579 TRUE
## 2 1 21.72480 2.251926 31.28994 296 4.247538 3.741706 TRUE
## 3 2 21.69482 2.190549 14.25503 79 3.396313 3.740829 FALSE
## 4 3 21.71821 2.153715 29.12643 178 3.439741 3.739949 FALSE
## 5 4 21.69484 2.116351 14.36555 274 3.463175 3.739067 FALSE
## 6 5 21.71804 2.078963 28.51150 321 3.267721 3.738181 FALSE
## 7 6 21.69647 2.046559 15.08181 260 3.232088 3.737291 FALSE
##
## attr(,"class")
## [1] "gofOutlier"lag_dat["241", "Glu_Gln_conc_molal"] <- NA
lag_dat["296", "Glu_Gln_conc_molal"] <- NA
#exclude PS17_047 from all analyses due to outlier status on multiple metabolites and issues in spectrum (noisy, structure in resid)
lag_dat["241", "NAA_conc_molal"] <- NA
lag_dat["241", "Cre_PCr_conc_molal"] <- NA
lag_dat["241", "Cho_GPC_PCh_conc_molal"] <- NA
write.csv(lag_dat, "/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/lag_data_molal_outliers_removed_n320.csv", row.names = FALSE)##Descriptive stats
#ACG - count individual subject N and N of females
acg_sub_sex <- acg_dat %>%
dplyr::select(subj_id, female)
acg_unique_subs<-unique(acg_sub_sex)
length(acg_unique_subs$subj_id)## [1] 102sum(acg_unique_subs$female)## [1] 47describe(acg_dat)## vars n mean sd median trimmed mad min
## study_code* 1 113 57.00 32.76 57.00 57.00 41.51 1.00
## subj_id* 2 113 51.42 28.89 54.00 51.46 35.58 1.00
## spec_id 3 113 5072.42 3034.48 3244.00 4708.09 1254.28 2058.00
## date_scan* 4 113 49.21 28.71 50.00 49.32 37.06 1.00
## hand* 5 113 4.54 1.04 5.00 4.78 0.00 1.00
## female 6 113 0.47 0.50 0.00 0.46 0.00 0.00
## mri_age_y 7 113 4.06 1.10 3.68 3.97 1.00 2.34
## folder* 8 113 57.00 32.76 57.00 57.00 41.51 1.00
## File_ID* 9 113 57.00 32.76 57.00 57.00 41.51 1.00
## fGM 10 113 0.77 0.05 0.78 0.78 0.05 0.55
## fWM 11 113 0.04 0.03 0.03 0.03 0.02 0.00
## fCSF 12 113 0.19 0.05 0.18 0.19 0.04 0.09
## NAA_conc_molal 13 113 14.89 1.32 15.03 14.94 1.25 11.24
## Cre_PCr_conc_molal 14 113 12.04 1.07 12.04 12.06 0.84 8.41
## Cho_GPC_PCh_conc_molal 15 112 2.50 0.34 2.48 2.49 0.30 1.57
## Ins_conc_molal 16 112 7.64 1.53 7.71 7.72 1.46 2.99
## Glu_conc_molal 17 113 22.38 2.25 22.50 22.46 2.01 16.03
## Glu_Gln_conc_molal 18 113 26.83 2.62 26.72 26.81 2.64 20.84
## max range skew kurtosis se
## study_code* 113.00 112.00 0.00 -1.23 3.08
## subj_id* 102.00 101.00 -0.05 -1.15 2.72
## spec_id 11448.00 9390.00 0.87 -0.81 285.46
## date_scan* 97.00 96.00 -0.04 -1.29 2.70
## hand* 6.00 5.00 -2.07 3.38 0.10
## female 1.00 1.00 0.12 -2.00 0.05
## mri_age_y 7.36 5.03 0.72 -0.49 0.10
## folder* 113.00 112.00 0.00 -1.23 3.08
## File_ID* 113.00 112.00 0.00 -1.23 3.08
## fGM 0.87 0.32 -1.14 1.59 0.01
## fWM 0.29 0.29 4.07 25.28 0.00
## fCSF 0.36 0.27 0.90 0.92 0.00
## NAA_conc_molal 18.30 7.05 -0.37 0.21 0.12
## Cre_PCr_conc_molal 14.76 6.35 -0.36 0.83 0.10
## Cho_GPC_PCh_conc_molal 3.33 1.76 0.16 0.08 0.03
## Ins_conc_molal 10.89 7.90 -0.51 0.35 0.14
## Glu_conc_molal 28.20 12.17 -0.31 0.40 0.21
## Glu_Gln_conc_molal 34.26 13.42 0.13 -0.01 0.25#LAG - count individual subject N and N of females
lag_sub_sex <- lag_dat %>%
dplyr::select(subj_id, female)
lag_unique_subs<-unique(lag_sub_sex)
length(lag_unique_subs$subj_id)## [1] 95sum(lag_unique_subs$female)## [1] 47describe(lag_dat)## vars n mean sd median trimmed mad min
## study_code* 1 321 160.94 92.72 161.00 161.00 118.61 1.00
## subj_id* 2 321 40.27 24.55 41.00 39.16 28.17 1.00
## spec_id 3 320 7744.26 3536.93 6379.50 7565.04 3846.61 3118.00
## date_scan* 4 321 143.73 82.33 143.00 143.92 103.78 1.00
## hand* 5 321 4.64 0.90 5.00 4.88 0.00 1.00
## female 6 321 0.47 0.50 0.00 0.46 0.00 0.00
## mri_age_y 7 321 5.95 1.90 5.39 5.79 1.79 2.41
## folder* 8 321 161.00 92.81 161.00 161.00 118.61 1.00
## File_ID* 9 321 161.00 92.81 161.00 161.00 118.61 1.00
## fGM 10 321 0.62 0.10 0.64 0.63 0.10 0.30
## fWM 11 321 0.35 0.11 0.34 0.34 0.11 0.00
## fCSF 12 321 0.03 0.02 0.03 0.03 0.01 0.00
## NAA_conc_molal 13 318 14.39 1.02 14.36 14.37 0.88 11.38
## Cre_PCr_conc_molal 14 320 10.25 0.91 10.17 10.23 0.83 7.28
## Cho_GPC_PCh_conc_molal 15 320 2.21 0.27 2.20 2.20 0.24 1.40
## Ins_conc_molal 16 319 6.37 1.19 6.37 6.39 1.08 2.04
## Glu_conc_molal 17 321 17.88 1.90 17.87 17.88 1.73 10.63
## Glu_Gln_conc_molal 18 319 21.69 2.19 21.74 21.72 2.02 14.26
## max range skew kurtosis se
## study_code* 320.00 319.00 0.00 -1.21 5.17
## subj_id* 95.00 94.00 0.28 -0.83 1.37
## spec_id 14205.00 11087.00 0.35 -1.37 197.72
## date_scan* 284.00 283.00 0.00 -1.19 4.60
## hand* 6.00 5.00 -2.35 4.48 0.05
## female 1.00 1.00 0.13 -1.99 0.03
## mri_age_y 11.13 8.72 0.66 -0.46 0.11
## folder* 321.00 320.00 0.00 -1.21 5.18
## File_ID* 321.00 320.00 0.00 -1.21 5.18
## fGM 0.80 0.50 -0.75 0.40 0.01
## fWM 0.66 0.65 0.43 0.20 0.01
## fCSF 0.21 0.21 3.12 16.87 0.00
## NAA_conc_molal 18.06 6.68 0.20 0.95 0.06
## Cre_PCr_conc_molal 14.09 6.81 0.17 1.33 0.05
## Cho_GPC_PCh_conc_molal 3.21 1.81 0.22 0.80 0.02
## Ins_conc_molal 11.78 9.74 -0.13 1.69 0.07
## Glu_conc_molal 24.87 14.24 -0.07 1.22 0.11
## Glu_Gln_conc_molal 29.13 14.87 -0.12 0.74 0.12#Calculate total Ns from both datasets combined
all_dat<-bind_rows(acg_dat, lag_dat)
length(unique(all_dat$File_ID))## [1] 434length(unique(all_dat$subj_id))## [1] 127all_sub_sex <- all_dat %>%
dplyr::select(subj_id, female)
all_unique_subs<-unique(all_sub_sex)
sum(all_unique_subs$female)## [1] 62describe(all_dat)## vars n mean sd median trimmed mad min
## study_code* 1 434 216.94 124.93 216.50 216.97 160.12 1.00
## subj_id* 2 434 53.87 33.15 54.00 52.20 38.55 1.00
## spec_id 3 433 7046.99 3606.23 5832.00 6818.88 3964.47 2058.00
## date_scan* 4 434 184.27 105.60 183.00 184.45 133.43 1.00
## hand* 5 434 4.62 0.93 5.00 4.86 0.00 1.00
## female 6 434 0.47 0.50 0.00 0.46 0.00 0.00
## mri_age_y 7 434 5.46 1.92 4.95 5.28 1.68 2.34
## folder* 8 434 217.04 125.00 217.50 217.05 160.12 1.00
## File_ID* 9 434 217.50 125.43 217.50 217.50 160.86 1.00
## fGM 10 434 0.66 0.11 0.67 0.67 0.11 0.30
## fWM 11 434 0.27 0.17 0.29 0.27 0.16 0.00
## fCSF 12 434 0.07 0.08 0.03 0.06 0.03 0.00
## NAA_conc_molal 13 431 14.52 1.13 14.48 14.51 0.99 11.24
## Cre_PCr_conc_molal 14 433 10.71 1.24 10.54 10.65 1.10 7.28
## Cho_GPC_PCh_conc_molal 15 432 2.28 0.32 2.26 2.27 0.26 1.40
## Ins_conc_molal 16 431 6.70 1.40 6.66 6.68 1.28 2.04
## Glu_conc_molal 17 434 19.05 2.81 18.53 18.89 2.53 10.63
## Glu_Gln_conc_molal 18 432 23.04 3.23 22.54 22.85 2.77 14.26
## max range skew kurtosis se
## study_code* 432.00 431.00 0.00 -1.21 6.00
## subj_id* 127.00 126.00 0.33 -0.78 1.59
## spec_id 14205.00 12147.00 0.45 -1.19 173.30
## date_scan* 366.00 365.00 0.00 -1.18 5.07
## hand* 6.00 5.00 -2.29 4.27 0.04
## female 1.00 1.00 0.13 -1.99 0.02
## mri_age_y 11.13 8.79 0.79 -0.09 0.09
## folder* 433.00 432.00 0.00 -1.20 6.00
## File_ID* 434.00 433.00 0.00 -1.21 6.02
## fGM 0.87 0.57 -0.60 0.10 0.01
## fWM 0.66 0.66 -0.09 -0.87 0.01
## fCSF 0.36 0.36 1.37 0.83 0.00
## NAA_conc_molal 18.30 7.05 0.10 0.57 0.05
## Cre_PCr_conc_molal 14.76 7.48 0.44 0.20 0.06
## Cho_GPC_PCh_conc_molal 3.33 1.93 0.45 0.70 0.02
## Ins_conc_molal 11.78 9.74 0.08 0.70 0.07
## Glu_conc_molal 28.20 17.57 0.49 0.13 0.13
## Glu_Gln_conc_molal 34.26 20.00 0.54 0.34 0.16LME modelling ACG
run mixed-effects models for ACG with participant as random effect, age, sex and age x sex interaction as fixed effects NAA run on acg_data_molal_outliers_removed_n113.csv
acg_naa_age<- lmer(NAA_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_naa_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 383.3 394.2 -187.6 375.3 109
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.53495 -0.39238 0.05273 0.38044 1.68117
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 1.1187 1.0577
## Residual 0.6121 0.7824
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.8379 0.4588 108.1511 30.16 <2e-16 ***
## mri_age_y 0.2513 0.1079 105.1446 2.33 0.0217 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.959acg_naa_age_sex<- lmer(NAA_conc_molal ~ mri_age_y + female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_naa_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 385.2 398.8 -187.6 375.2 108
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.54531 -0.40334 0.05382 0.37656 1.67013
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 1.1189 1.0578
## Residual 0.6113 0.7819
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.85437 0.46506 109.75797 29.791 <2e-16 ***
## mri_age_y 0.25366 0.10837 104.64089 2.341 0.0211 *
## female -0.05648 0.26075 99.22276 -0.217 0.8289
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.926
## female -0.165 -0.098anova(acg_naa_age,acg_naa_age_sex)## Data: acg_dat
## Models:
## acg_naa_age: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_naa_age_sex: NAA_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_naa_age 4 383.25 394.16 -187.62 375.25
## acg_naa_age_sex 5 385.20 398.84 -187.60 375.20 0.0469 1 0.8285acg_naa_ageXsex<- lmer(NAA_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_naa_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 |
## subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 387.1 403.4 -187.5 375.1 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.55254 -0.39337 0.03275 0.37793 1.69016
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 1.1179 1.0573
## Residual 0.6105 0.7813
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.6886 0.6657 106.4622 20.562 <2e-16 ***
## mri_age_y 0.2954 0.1616 107.8036 1.827 0.0704 .
## female 0.2522 0.9254 110.1751 0.273 0.7857
## mri_age_y:female -0.0757 0.2178 109.0452 -0.348 0.7288
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.964
## female -0.719 0.694
## mri_g_y:fml 0.716 -0.742 -0.960anova(acg_naa_age,acg_naa_ageXsex)## Data: acg_dat
## Models:
## acg_naa_age: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_naa_ageXsex: NAA_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_naa_age 4 383.25 394.16 -187.62 375.25
## acg_naa_ageXsex 6 387.08 403.45 -187.54 375.08 0.1677 2 0.9196#including sex in model does not significantly improve fit Cre
acg_cre_age<- lmer(Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_cre_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 340.9 351.8 -166.5 332.9 109
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.12075 -0.39494 0.02556 0.37255 1.62021
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.6867 0.8287
## Residual 0.4877 0.6983
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 11.67366 0.38192 110.29328 30.566 <2e-16 ***
## mri_age_y 0.09221 0.08992 107.78610 1.026 0.307
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.960acg_cre_age_sex<- lmer(Cre_PCr_conc_molal ~ mri_age_y + female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_cre_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 342.9 356.5 -166.4 332.9 108
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1065 -0.3932 0.0080 0.3605 1.6020
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.6866 0.8286
## Residual 0.4871 0.6979
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 11.65746 0.38678 111.48623 30.139 <2e-16 ***
## mri_age_y 0.08988 0.09034 107.34558 0.995 0.322
## female 0.05586 0.21445 96.67410 0.260 0.795
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.927
## female -0.160 -0.100anova(acg_cre_age,acg_cre_age_sex)## Data: acg_dat
## Models:
## acg_cre_age: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_cre_age_sex: Cre_PCr_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_cre_age 4 340.94 351.85 -166.47 332.94
## acg_cre_age_sex 5 342.87 356.51 -166.44 332.87 0.0678 1 0.7945acg_cre_ageXsex<- lmer(Cre_PCr_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_cre_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + female + mri_age_y:female +
## (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 343.6 360.0 -165.8 331.6 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.11740 -0.37592 0.04868 0.38312 1.59263
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.7009 0.8372
## Residual 0.4638 0.6810
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 12.09093 0.54699 105.10134 22.104 <2e-16 ***
## mri_age_y -0.01931 0.13288 106.80701 -0.145 0.885
## female -0.77294 0.76525 111.06817 -1.010 0.315
## mri_age_y:female 0.20325 0.18016 110.02199 1.128 0.262
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.965
## female -0.715 0.689
## mri_g_y:fml 0.711 -0.738 -0.960anova(acg_cre_age,acg_cre_ageXsex)## Data: acg_dat
## Models:
## acg_cre_age: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_cre_ageXsex: Cre_PCr_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_cre_age 4 340.94 351.85 -166.47 332.94
## acg_cre_ageXsex 6 343.61 359.98 -165.81 331.61 1.3273 2 0.515Cho
acg_cho_age<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_cho_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 78.1 89.0 -35.0 70.1 108
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1325 -0.5979 -0.1286 0.6106 2.4400
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.03204 0.1790
## Residual 0.07887 0.2808
## Number of obs: 112, groups: subj_id, 101
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.74592 0.12031 111.92574 22.82 <2e-16 ***
## mri_age_y -0.06423 0.02842 111.94035 -2.26 0.0258 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.963acg_cho_age_sex<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_cho_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 78.2 91.8 -34.1 68.2 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.17836 -0.63945 -0.09041 0.54644 2.50032
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.03214 0.1793
## Residual 0.07694 0.2774
## Number of obs: 112, groups: subj_id, 101
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.77074 0.12063 111.75493 22.970 <2e-16 ***
## mri_age_y -0.06004 0.02834 111.86929 -2.118 0.0364 *
## female -0.09017 0.06493 103.29750 -1.389 0.1679
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.931
## female -0.148 -0.107anova(acg_cho_age,acg_cho_age_sex)## Data: acg_dat
## Models:
## acg_cho_age: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_cho_age_sex: Cho_GPC_PCh_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_cho_age 4 78.077 88.951 -35.039 70.077
## acg_cho_age_sex 5 78.166 91.759 -34.083 68.166 1.9111 1 0.1668acg_cho_ageXsex<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_cho_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + female + mri_age_y:female +
## (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 79.7 96.1 -33.9 67.7 106
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1308 -0.6109 -0.1102 0.5606 2.5638
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.03297 0.1816
## Residual 0.07579 0.2753
## Number of obs: 112, groups: subj_id, 101
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.69546 0.16726 107.68134 16.116 <2e-16 ***
## mri_age_y -0.04107 0.04072 109.53005 -1.008 0.315
## female 0.05961 0.23991 111.98464 0.248 0.804
## mri_age_y:female -0.03674 0.05661 111.93901 -0.649 0.518
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.965
## female -0.697 0.673
## mri_g_y:fml 0.694 -0.719 -0.963anova(acg_cho_age,acg_cho_ageXsex)## Data: acg_dat
## Models:
## acg_cho_age: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_cho_ageXsex: Cho_GPC_PCh_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_cho_age 4 78.077 88.951 -35.039 70.077
## acg_cho_ageXsex 6 79.747 96.058 -33.874 67.747 2.3297 2 0.312#including sex in model does not significantly improve fit Ins
acg_ins_age<- lmer(Ins_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_ins_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 416.5 427.3 -204.2 408.5 108
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.51194 -0.47558 0.09217 0.51196 1.57868
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.9268 0.9627
## Residual 1.3749 1.1726
## Number of obs: 112, groups: subj_id, 101
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.4742 0.5464 111.9652 15.509 <2e-16 ***
## mri_age_y -0.2049 0.1294 111.4483 -1.583 0.116
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.962acg_ins_age_sex<- lmer(Ins_conc_molal ~ mri_age_y + female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_ins_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 417.1 430.7 -203.5 407.1 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.47103 -0.48748 0.05729 0.54401 1.51171
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.8891 0.9429
## Residual 1.3819 1.1755
## Number of obs: 112, groups: subj_id, 101
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.5595 0.5485 111.9724 15.605 <2e-16 ***
## mri_age_y -0.1859 0.1296 111.3196 -1.435 0.154
## female -0.3489 0.2981 99.5526 -1.171 0.245
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.930
## female -0.140 -0.118anova(acg_ins_age,acg_ins_age_sex)## Data: acg_dat
## Models:
## acg_ins_age: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_ins_age_sex: Ins_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_ins_age 4 416.45 427.33 -204.23 408.45
## acg_ins_age_sex 5 417.09 430.69 -203.55 407.09 1.3598 1 0.2436acg_ins_ageXsex<- lmer(Ins_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_ins_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 |
## subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 418.2 434.5 -203.1 406.2 106
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.55603 -0.53421 0.07378 0.57588 1.65390
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.7615 0.8727
## Residual 1.4790 1.2162
## Number of obs: 112, groups: subj_id, 101
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 9.0615 0.7667 105.7180 11.819 <2e-16 ***
## mri_age_y -0.3138 0.1882 108.1469 -1.667 0.0983 .
## female -1.3693 1.0920 111.9981 -1.254 0.2125
## mri_age_y:female 0.2511 0.2588 111.8723 0.970 0.3340
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.965
## female -0.702 0.678
## mri_g_y:fml 0.702 -0.727 -0.963anova(acg_ins_age,acg_ins_ageXsex)## Data: acg_dat
## Models:
## acg_ins_age: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_ins_ageXsex: Ins_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_ins_age 4 416.45 427.33 -204.23 408.45
## acg_ins_ageXsex 6 418.20 434.51 -203.10 406.20 2.2549 2 0.3239Glx
acg_glx_age<- lmer(Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_glx_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 541.1 552.0 -266.5 533.1 109
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.48535 -0.37243 -0.00651 0.43846 1.68169
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 4.706 2.169
## Residual 2.330 1.526
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 25.7502 0.9206 106.5717 27.972 <2e-16 ***
## mri_age_y 0.2612 0.2163 103.0785 1.207 0.23
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.959acg_glx_age_sex<- lmer(Glu_Gln_conc_molal ~ mri_age_y + female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_glx_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 543.0 556.6 -266.5 533.0 108
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.47879 -0.39442 -0.02857 0.42851 1.70104
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 4.692 2.166
## Residual 2.335 1.528
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 25.7002 0.9334 108.5734 27.535 <2e-16 ***
## mri_age_y 0.2542 0.2173 102.5984 1.170 0.245
## female 0.1700 0.5257 98.8828 0.323 0.747
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.925
## female -0.167 -0.097anova(acg_glx_age,acg_glx_age_sex)## Data: acg_dat
## Models:
## acg_glx_age: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_glx_age_sex: Glu_Gln_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_glx_age 4 541.09 552.00 -266.55 533.09
## acg_glx_age_sex 5 542.99 556.62 -266.49 532.99 0.1045 1 0.7465acg_glx_ageXsex<- lmer(Glu_Gln_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=acg_dat, REML=FALSE)
summary(acg_glx_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + female + mri_age_y:female +
## (1 | subj_id)
## Data: acg_dat
##
## AIC BIC logLik deviance df.resid
## 542.5 558.9 -265.2 530.5 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.53333 -0.43763 0.07324 0.38166 1.63409
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 4.573 2.139
## Residual 2.297 1.515
## Number of obs: 113, groups: subj_id, 102
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 27.2150 1.3265 106.5825 20.516 <2e-16 ***
## mri_age_y -0.1271 0.3220 107.8777 -0.395 0.694
## female -2.6328 1.8384 109.4802 -1.432 0.155
## mri_age_y:female 0.6875 0.4325 108.2254 1.590 0.115
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.964
## female -0.722 0.696
## mri_g_y:fml 0.718 -0.745 -0.959anova(acg_glx_age,acg_glx_ageXsex)## Data: acg_dat
## Models:
## acg_glx_age: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_glx_ageXsex: Glu_Gln_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_glx_age 4 541.09 552.00 -266.55 533.09
## acg_glx_ageXsex 6 542.49 558.85 -265.25 530.49 2.603 2 0.2721Plot results ACG
Create scatter plot with regression line based on model predictions NAA
acg_dat <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/acg_data_molal_outliers_removed_n113.csv")
#pull intercept and slope from model for regression line
fixef.acg_naa<-fixef(acg_naa_age)
#plot
plot_acg_naa_age<-ggplot(acg_dat, aes(x=mri_age_y, y=NAA_conc_molal,color=as.factor(female)))
plot_acg_naa_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(title = "Anterior Cingulate", x = "Age (years)", y = "tNAA (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.acg_naa[[1]], slope= fixef.acg_naa[[2]])) +
theme_classic()
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/acg_naa_age_viridis.png")## Saving 7 x 5 in imageCho
#pull intercept and slope from model for regression line
fixef.acg_cho<-fixef(acg_cho_age)
#plot
plot_acg_cho_age<-ggplot(acg_dat, aes(x=mri_age_y, y=Cho_GPC_PCh_conc_molal,color=as.factor(female)))
plot_acg_cho_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(title = "Anterior Cingulate", x = "Age (years)", y = "tCho (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.acg_cho[[1]], slope= fixef.acg_cho[[2]])) +
theme_classic()## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 row(s) containing missing values (geom_path).
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/acg_cho_age_viridis.png")## Saving 7 x 5 in image## Warning: Removed 1 rows containing missing values (geom_point).
## Warning: Removed 1 row(s) containing missing values (geom_path).Ins
#pull intercept and slope from model for regression line
fixef.acg_ins<-fixef(acg_ins_age)
#plot
plot_acg_ins_age<-ggplot(acg_dat, aes(x=mri_age_y, y=Ins_conc_molal,color=as.factor(female)))
plot_acg_ins_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(title = "Anterior Cingulate", x = "Age (years)", y = "Ins (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.acg_ins[[1]], slope= fixef.acg_ins[[2]])) +
theme_classic()## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 row(s) containing missing values (geom_path).
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/acg_ins_age_viridis.png")## Saving 7 x 5 in image## Warning: Removed 1 rows containing missing values (geom_point).
## Warning: Removed 1 row(s) containing missing values (geom_path).LME modelling LAG
run mixed-effects models for lag with participant as random effect, age, sex and age x sex interaction as fixed effects NAA, substantial longitudinal data so we include random slope run on lag_data_molal_outliers_removed_n320.csv
#basic model effect of age on intercept
lag_naa_age<- lmer(NAA_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_naa_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 896 911 -444 888 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4313 -0.5792 -0.0185 0.5440 3.6263
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1641 0.4051
## Residual 0.8272 0.9095
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.82375 0.18588 309.72511 74.369 < 2e-16 ***
## mri_age_y 0.09699 0.02902 316.51634 3.343 0.00093 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.929lag_naa_age_sex<- lmer(NAA_conc_molal ~ mri_age_y + female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_naa_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 898.0 916.8 -444.0 888.0 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4286 -0.5791 -0.0185 0.5452 3.6231
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1637 0.4046
## Residual 0.8275 0.9097
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.827819 0.194854 290.918032 70.965 < 2e-16 ***
## mri_age_y 0.097086 0.029038 316.589446 3.343 0.000927 ***
## female -0.009878 0.137557 93.706753 -0.072 0.942903
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.874
## female -0.300 -0.039anova(lag_naa_age,lag_naa_age_sex)## Data: lag_dat
## Models:
## lag_naa_age: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_naa_age_sex: NAA_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_naa_age 4 895.96 911.01 -443.98 887.96
## lag_naa_age_sex 5 897.95 916.76 -443.98 887.95 0.0051 1 0.943lag_naa_ageXsex<- lmer(NAA_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_naa_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 |
## subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 899.9 922.5 -444.0 887.9 312
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4372 -0.5859 -0.0225 0.5510 3.6185
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1645 0.4056
## Residual 0.8268 0.9093
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.79473 0.24812 314.11810 55.598 < 2e-16 ***
## mri_age_y 0.10272 0.03910 313.83715 2.627 0.00903 **
## female 0.06577 0.37459 308.18552 0.176 0.86075
## mri_age_y:female -0.01265 0.05838 316.98288 -0.217 0.82864
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.924
## female -0.662 0.612
## mri_g_y:fml 0.619 -0.670 -0.930anova(lag_naa_age,lag_naa_ageXsex)## Data: lag_dat
## Models:
## lag_naa_age: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_naa_ageXsex: NAA_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_naa_age 4 895.96 911.01 -443.98 887.96
## lag_naa_ageXsex 6 899.91 922.48 -443.95 887.91 0.0518 2 0.9744#including sex in model does not significantly improve fit
#test quadratic effect of age
#create age^2 variable
lag_dat$mri_age_y2 <- (lag_dat$mri_age_y)^2
lag_naa_age2<- lmer(NAA_conc_molal ~ mri_age_y + mri_age_y2 + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_naa_age2) #age^2 p=.0447## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + mri_age_y2 + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 893.9 912.7 -442.0 883.9 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5451 -0.5819 -0.0270 0.4941 3.5330
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1700 0.4123
## Residual 0.8119 0.9011
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 12.70181 0.58588 310.95949 21.680 <2e-16 ***
## mri_age_y 0.47130 0.18749 303.10334 2.514 0.0125 *
## mri_age_y2 -0.02836 0.01402 297.51463 -2.023 0.0440 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.983
## mri_age_y2 0.949 -0.988anova(lag_naa_age, lag_naa_age2) #fit improved p=.044## Data: lag_dat
## Models:
## lag_naa_age: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_naa_age2: NAA_conc_molal ~ mri_age_y + mri_age_y2 + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_naa_age 4 895.96 911.01 -443.98 887.96
## lag_naa_age2 5 893.91 912.72 -441.96 883.91 4.0455 1 0.04429 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Cre
lag_cre_age<- lmer(Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cre_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 848.2 863.3 -420.1 840.2 316
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9764 -0.5307 -0.1261 0.5525 4.1548
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.08477 0.2912
## Residual 0.73683 0.8584
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 10.49912 0.16805 308.49608 62.476 <2e-16 ***
## mri_age_y -0.04183 0.02646 319.87112 -1.581 0.115
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.937lag_cre_age_sex<- lmer(Cre_PCr_conc_molal ~ mri_age_y + female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cre_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 849.7 868.6 -419.9 839.7 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9528 -0.5658 -0.1208 0.5434 4.1925
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.08478 0.2912
## Residual 0.73557 0.8577
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 10.46482 0.17478 289.38001 59.873 <2e-16 ***
## mri_age_y -0.04267 0.02646 319.88902 -1.612 0.108
## female 0.08318 0.11739 84.73092 0.709 0.481
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.887
## female -0.277 -0.044anova(lag_cre_age,lag_cre_age_sex)## Data: lag_dat
## Models:
## lag_cre_age: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_cre_age_sex: Cre_PCr_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_cre_age 4 848.23 863.31 -420.12 840.23
## lag_cre_age_sex 5 849.73 868.57 -419.87 839.73 0.5017 1 0.4788lag_cre_ageXsex<- lmer(Cre_PCr_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cre_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + female + mri_age_y:female +
## (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 850.0 872.6 -419.0 838.0 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9579 -0.5653 -0.0910 0.5718 4.2323
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.08177 0.2860
## Residual 0.73339 0.8564
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 10.65304 0.22484 314.11103 47.381 <2e-16 ***
## mri_age_y -0.07476 0.03585 318.90184 -2.086 0.0378 *
## female -0.33761 0.33688 307.33953 -1.002 0.3171
## mri_age_y:female 0.07047 0.05297 319.97938 1.330 0.1844
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.935
## female -0.667 0.624
## mri_g_y:fml 0.632 -0.677 -0.938anova(lag_cre_age,lag_cre_ageXsex)## Data: lag_dat
## Models:
## lag_cre_age: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_cre_ageXsex: Cre_PCr_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_cre_age 4 848.23 863.31 -420.12 840.23
## lag_cre_ageXsex 6 849.97 872.58 -418.98 837.97 2.263 2 0.3226Cho
lag_cho_age<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 51.3 66.4 -21.7 43.3 316
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6539 -0.5810 -0.0562 0.4659 4.1202
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.0149 0.1221
## Residual 0.0561 0.2368
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.346670 0.049038 310.144986 47.854 < 2e-16 ***
## mri_age_y -0.023315 0.007581 315.773363 -3.075 0.00229 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.921lag_cho_age_sex<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 51.6 70.5 -20.8 41.6 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6317 -0.5893 -0.0587 0.4688 4.1717
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.01469 0.1212
## Residual 0.05586 0.2364
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.36786 0.05147 283.09492 46.006 < 2e-16 ***
## mri_age_y -0.02286 0.00757 315.96555 -3.019 0.00274 **
## female -0.05024 0.03808 80.73098 -1.319 0.19083
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.860
## female -0.312 -0.046anova(lag_cho_age,lag_cho_age_sex)## Data: lag_dat
## Models:
## lag_cho_age: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_cho_age_sex: Cho_GPC_PCh_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_cho_age 4 51.345 66.419 -21.673 43.345
## lag_cho_age_sex 5 51.610 70.452 -20.805 41.610 1.735 1 0.1878lag_cho_ageXsex<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + female + mri_age_y:female +
## (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 49.9 72.5 -18.9 37.9 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6639 -0.5411 -0.0413 0.4642 4.0526
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.01452 0.1205
## Residual 0.05521 0.2350
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.28998 0.06505 314.68723 35.206 <2e-16 ***
## mri_age_y -0.00954 0.01019 312.12139 -0.936 0.3497
## female 0.12488 0.09791 309.31059 1.275 0.2031
## mri_age_y:female -0.02931 0.01511 316.65469 -1.939 0.0533 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.916
## female -0.664 0.609
## mri_g_y:fml 0.617 -0.674 -0.922anova(lag_cho_age,lag_cho_ageXsex)## Data: lag_dat
## Models:
## lag_cho_age: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_cho_ageXsex: Cho_GPC_PCh_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_cho_age 4 51.345 66.419 -21.673 43.345
## lag_cho_ageXsex 6 49.871 72.481 -18.936 37.871 5.4741 2 0.06476 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#including sex in model does not significantly improve fit Ins
lag_ins_age<- lmer(Ins_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_ins_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1013.3 1028.3 -502.6 1005.3 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3984 -0.5594 -0.0111 0.5911 4.6516
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1768 0.4205
## Residual 1.2227 1.1057
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.28473 0.21938 307.28604 28.648 <2e-16 ***
## mri_age_y 0.01632 0.03444 318.45543 0.474 0.636
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.934lag_ins_age_sex<- lmer(Ins_conc_molal ~ mri_age_y + female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_ins_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1015.3 1034.1 -502.6 1005.3 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.4003 -0.5588 -0.0120 0.5901 4.6506
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1765 0.4201
## Residual 1.2229 1.1058
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.286355 0.228788 284.644902 27.477 <2e-16 ***
## mri_age_y 0.016356 0.034473 318.479871 0.474 0.635
## female -0.003854 0.157311 78.609245 -0.024 0.981
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.881
## female -0.284 -0.046anova(lag_ins_age,lag_ins_age_sex)## Data: lag_dat
## Models:
## lag_ins_age: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_ins_age_sex: Ins_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_ins_age 4 1013.3 1028.3 -502.64 1005.3
## lag_ins_age_sex 5 1015.3 1034.1 -502.64 1005.3 6e-04 1 0.9807lag_ins_ageXsex<- lmer(Ins_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_ins_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 |
## subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1017.1 1039.7 -502.6 1005.1 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.4248 -0.5418 -0.0121 0.5899 4.6348
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1767 0.4204
## Residual 1.2221 1.1055
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.361751 0.294417 312.731852 21.608 <2e-16 ***
## mri_age_y 0.003458 0.046822 316.620308 0.074 0.941
## female -0.171702 0.441468 305.644875 -0.389 0.698
## mri_age_y:female 0.028161 0.069170 318.668940 0.407 0.684
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.930
## female -0.667 0.620
## mri_g_y:fml 0.630 -0.677 -0.934anova(lag_ins_age,lag_ins_ageXsex)## Data: lag_dat
## Models:
## lag_ins_age: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_ins_ageXsex: Ins_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_ins_age 4 1013.3 1028.3 -502.64 1005.3
## lag_ins_ageXsex 6 1017.1 1039.7 -502.56 1005.1 0.1663 2 0.9202Glx
lag_glx_age<- lmer(Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularisSingular(lag_glx_age) # TRUE: model is almost/near singular (parameters are on the boundary of the feasible parameter space: random effects variance = 0)## [1] TRUEsummary(lag_glx_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1396.1 1411.1 -694.0 1388.1 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6700 -0.6848 0.0127 0.6729 3.5338
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.000 0.000
## Residual 4.542 2.131
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 23.23676 0.39326 319.00000 59.088 < 2e-16 ***
## mri_age_y -0.25961 0.06309 319.00000 -4.115 4.93e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.953
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingularlag_glx_age_sex<- lmer(Glu_Gln_conc_molal ~ mri_age_y + female + (1 | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_glx_age_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + female + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1396.7 1415.5 -693.4 1386.7 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7491 -0.6907 0.0192 0.7055 3.6028
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.000 0.000
## Residual 4.523 2.127
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 23.12126 0.40478 319.00000 57.120 < 2e-16 ***
## mri_age_y -0.26202 0.06299 319.00000 -4.160 4.1e-05 ***
## female 0.27789 0.23880 319.00000 1.164 0.245
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.915
## female -0.245 -0.033
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingularanova(lag_glx_age,lag_glx_age_sex)## Data: lag_dat
## Models:
## lag_glx_age: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_glx_age_sex: Glu_Gln_conc_molal ~ mri_age_y + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_glx_age 4 1396.1 1411.1 -694.04 1388.1
## lag_glx_age_sex 5 1396.7 1415.5 -693.36 1386.7 1.3513 1 0.245lag_glx_ageXsex<- lmer(Glu_Gln_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_glx_ageXsex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + female + mri_age_y:female +
## (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1397.8 1420.4 -692.9 1385.8 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7059 -0.6942 -0.0013 0.6785 3.6197
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.00 0.000
## Residual 4.51 2.124
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 23.44862 0.52588 319.00000 44.589 < 2e-16 ***
## mri_age_y -0.31768 0.08502 319.00000 -3.737 0.000221 ***
## female -0.45362 0.78868 319.00000 -0.575 0.565583
## mri_age_y:female 0.12295 0.12636 319.00000 0.973 0.331265
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y female
## mri_age_y -0.951
## female -0.667 0.634
## mri_g_y:fml 0.640 -0.673 -0.953
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingularanova(lag_glx_age,lag_glx_ageXsex)## Data: lag_dat
## Models:
## lag_glx_age: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_glx_ageXsex: Glu_Gln_conc_molal ~ mri_age_y + female + mri_age_y:female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_glx_age 4 1396.1 1411.1 -694.04 1388.1
## lag_glx_ageXsex 6 1397.8 1420.4 -692.89 1385.8 2.2967 2 0.3172Plot results LAG
Create scatter plot with regression line based on model predictions NAA
lag_dat<-read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/lag_data_molal_outliers_removed_n320.csv")
#pull intercept and slope from model for regression line
fixef.lag_naa<-fixef(lag_naa_age)
#plot
plot_lag_naa_age<-ggplot(lag_dat, aes(x=mri_age_y, y=NAA_conc_molal,color=as.factor(female)))
plot_lag_naa_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(title = "Left Temporo-Parietal", x = "Age (years)", y = "tNAA (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.lag_naa[[1]], slope= fixef.lag_naa[[2]])) +
theme_classic()## Warning: Removed 3 rows containing missing values (geom_point).## Warning: Removed 3 row(s) containing missing values (geom_path).
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/lag_naa_age_viridis.png")## Saving 7 x 5 in image## Warning: Removed 3 rows containing missing values (geom_point).
## Warning: Removed 3 row(s) containing missing values (geom_path).Cho
#pull intercept and slope from model for regression line
fixef.lag_cho<-fixef(lag_cho_age)
#plot
plot_lag_cho_age<-ggplot(lag_dat, aes(x=mri_age_y, y=Cho_GPC_PCh_conc_molal,color=as.factor(female)))
plot_lag_cho_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(title = "Left Temporo-Parietal", x = "Age (years)", y = "tCho (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.lag_cho[[1]], slope= fixef.lag_cho[[2]])) +
theme_classic()## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 row(s) containing missing values (geom_path).
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/lag_cho_age_viridis.png")## Saving 7 x 5 in image## Warning: Removed 1 rows containing missing values (geom_point).
## Warning: Removed 1 row(s) containing missing values (geom_path).Glx
#pull intercept and slope from model for regression line
fixef.lag_glx<-fixef(lag_glx_age)
#plot
plot_lag_glx_age<-ggplot(lag_dat, aes(x=mri_age_y, y=Glu_Gln_conc_molal,color=as.factor(female)))
plot_lag_glx_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(title = "Left Temporo-Parietal", x = "Age (years)", y = "Glx (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.lag_glx[[1]], slope= fixef.lag_glx[[2]])) +
theme_classic()## Warning: Removed 2 rows containing missing values (geom_point).## Warning: Removed 2 row(s) containing missing values (geom_path).
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/lag_glx_age_viridis.png")## Saving 7 x 5 in image## Warning: Removed 2 rows containing missing values (geom_point).
## Warning: Removed 2 row(s) containing missing values (geom_path).##Combine results plots into 1 figure Simplify plots to legend & region don’t repeat in each panel
library(patchwork)
p1_acg_naa_age<-plot_acg_naa_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(title = "Anterior Cingulate", x = "Age (years)", y = "tNAA (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.acg_naa[[1]], slope= fixef.acg_naa[[2]])) +
theme_classic() +
theme(plot.title= element_text(hjust = 0.5))
p4_lag_naa_age<-plot_lag_naa_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(title = "Left Temporo-Parietal", x = "Age (years)", y = "tNAA (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
theme_classic() +
#theme(legend.position = c(0.7, 0.95), legend.background = element_rect(fill='transparent')) +
geom_abline(aes(intercept = fixef.lag_naa[[1]], slope= fixef.lag_naa[[2]])) +
theme(plot.title= element_text(hjust = 0.5))
p2_acg_cho_age<-plot_acg_cho_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(x = "Age (years)", y = "tCho (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.acg_cho[[1]], slope= fixef.acg_cho[[2]])) +
theme_classic()
p5_lag_cho_age<-plot_lag_cho_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(x = "Age (years)", y = "tCho (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.lag_cho[[1]], slope= fixef.lag_cho[[2]])) +
theme_classic()
p3_acg_ins_age<-plot_acg_ins_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(x = "Age (years)", y = "Ins (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.acg_ins[[1]], slope= fixef.acg_ins[[2]])) +
theme_classic()
p6_lag_glx_age<-plot_lag_glx_age +
geom_point() +
geom_line(aes(group = as.factor(subj_id))) +
labs(x = "Age (years)", y = "Glx (molal)") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.lag_glx[[1]], slope= fixef.lag_glx[[2]])) +
theme_classic()
results_fig<-p1_acg_naa_age + p4_lag_naa_age + p2_acg_cho_age + p5_lag_cho_age + p3_acg_ins_age + p6_lag_glx_age + plot_layout(nrow=3, guides ='collect')
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/results_fig_grid_viridis.png", results_fig, width = 2000, units = "px")## Saving 2000 x 1500 px image## Warning: Removed 3 rows containing missing values (geom_point).## Warning: Removed 3 row(s) containing missing values (geom_path).## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 row(s) containing missing values (geom_path).## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 row(s) containing missing values (geom_path).## Warning: Removed 1 rows containing missing values (geom_point).## Warning: Removed 1 row(s) containing missing values (geom_path).## Warning: Removed 2 rows containing missing values (geom_point).## Warning: Removed 2 row(s) containing missing values (geom_path).##Correct for multiple comparisons Use FDR method (Benjamini & Hochberg) to correct for mutliple comparisons (5 metabolites x 2 voxels = 10 models)
#list models (voxel/metabolite pair)
models<- c("acg_naa", "acg_cr", "acg_cho", "acg_ins", "acg_glx", "lag_naa", "lag_cr", "lag_cho", "lag_ins", "lag_glx")
#list p-values for each model, matching order of models
p<-c(.0217, .307, .0258, .0494, .23, .00093, .115, .00229, .636, 0)
fdr_table<-data.frame(models, p)
fdr_table$p_fdr<-p.adjust(p, method = "fdr", n= length(p))
fdr_table## models p p_fdr
## 1 acg_naa 0.02170 0.051600000
## 2 acg_cr 0.30700 0.341111111
## 3 acg_cho 0.02580 0.051600000
## 4 acg_ins 0.04940 0.082333333
## 5 acg_glx 0.23000 0.287500000
## 6 lag_naa 0.00093 0.004650000
## 7 lag_cr 0.11500 0.164285714
## 8 lag_cho 0.00229 0.007633333
## 9 lag_ins 0.63600 0.636000000
## 10 lag_glx 0.00000 0.000000000#run FDR without LAG-Glx in case it's invalid:
p2<-c(.0217, .307, .0258, .0494, .23, .00093, .115, .00229, .636)
p.adjust(p2, method = "fdr", n= length(p2))## [1] 0.0580500 0.3453750 0.0580500 0.0889200 0.2957143 0.0083700 0.1725000
## [8] 0.0103050 0.6360000##Additional models Next steps following OHBM abstract submission (Dec. 2021): test random slopes (LAG only)
#read in data
lag_dat<-read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/lag_data_molal_outliers_removed_n320.csv")
#basic age model NAA
lag_naa_age<- lmer(NAA_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_naa_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 896 911 -444 888 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4313 -0.5792 -0.0185 0.5440 3.6263
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1641 0.4051
## Residual 0.8272 0.9095
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.82375 0.18588 309.72511 74.369 < 2e-16 ***
## mri_age_y 0.09699 0.02902 316.51634 3.343 0.00093 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.929#test random slope NAA
lag_naa_age_randslope <- lmer(NAA_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id), data=lag_dat, REML=FALSE)## Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
## Model failed to converge with max|grad| = 0.00298042 (tol = 0.002, component 1)summary(lag_naa_age_randslope)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 899.3 921.9 -443.7 887.3 312
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4425 -0.5829 -0.0345 0.5457 3.6775
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## subj_id (Intercept) 0.438346 0.66208
## mri_age_y 0.005125 0.07159 -0.81
## Residual 0.808385 0.89910
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.81838 0.19532 84.82239 70.749 < 2e-16 ***
## mri_age_y 0.09843 0.03011 76.83548 3.269 0.00162 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.936
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00298042 (tol = 0.002, component 1)#model failed to converge
#basic age model Cre
lag_cre_age<- lmer(Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cre_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 848.2 863.3 -420.1 840.2 316
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9764 -0.5307 -0.1261 0.5525 4.1548
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.08477 0.2912
## Residual 0.73683 0.8584
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 10.49912 0.16805 308.49608 62.476 <2e-16 ***
## mri_age_y -0.04183 0.02646 319.87112 -1.581 0.115
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.937#test random slope Cre
lag_cre_age_randslope<- lmer(Cre_PCr_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id), data=lag_dat, REML=FALSE)## Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
## Model failed to converge with max|grad| = 0.00521307 (tol = 0.002, component 1)summary(lag_cre_age_randslope)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 848.6 871.2 -418.3 836.6 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1240 -0.5211 -0.1430 0.5938 4.3226
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## subj_id (Intercept) 0.594776 0.77122
## mri_age_y 0.008986 0.09479 -0.96
## Residual 0.705980 0.84023
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 10.50247 0.18591 91.11070 56.492 <2e-16 ***
## mri_age_y -0.04111 0.02825 82.73282 -1.455 0.149
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.951
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00521307 (tol = 0.002, component 1)#model failed to converge
#basic age model Cho
lag_cho_age<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 51.3 66.4 -21.7 43.3 316
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6539 -0.5810 -0.0562 0.4659 4.1202
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.0149 0.1221
## Residual 0.0561 0.2368
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.346670 0.049038 310.144986 47.854 < 2e-16 ***
## mri_age_y -0.023315 0.007581 315.773363 -3.075 0.00229 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.921#test random slope Cho
lag_cho_age_randslope<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_age_randslope)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 55.0 77.6 -21.5 43.0 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6693 -0.5804 -0.0496 0.4832 4.0955
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## subj_id (Intercept) 0.0211331 0.14537
## mri_age_y 0.0002882 0.01698 -0.56
## Residual 0.0550144 0.23455
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.350155 0.050009 57.473557 46.994 < 2e-16 ***
## mri_age_y -0.024037 0.007888 61.859198 -3.047 0.00339 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.924#compare models
anova(lag_cho_age, lag_cho_age_randslope) #fit did not improve with inclusion of random slopes## Data: lag_dat
## Models:
## lag_cho_age: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_cho_age_randslope: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_cho_age 4 51.345 66.419 -21.673 43.345
## lag_cho_age_randslope 6 55.014 77.624 -21.507 43.014 0.3314 2 0.8473#basic age model Ins
lag_ins_age<- lmer(Ins_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_ins_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1013.3 1028.3 -502.6 1005.3 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3984 -0.5594 -0.0111 0.5911 4.6516
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1768 0.4205
## Residual 1.2227 1.1057
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.28473 0.21938 307.28604 28.648 <2e-16 ***
## mri_age_y 0.01632 0.03444 318.45543 0.474 0.636
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.934#test random slope Ins
lag_ins_age_randslope<- lmer(Ins_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_ins_age_randslope)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1016.6 1039.1 -502.3 1004.6 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2826 -0.5584 -0.0263 0.5849 4.6042
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## subj_id (Intercept) 0.414651 0.64393
## mri_age_y 0.001291 0.03593 -1.00
## Residual 1.211663 1.10076
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.29900 0.22582 145.84527 27.894 <2e-16 ***
## mri_age_y 0.01355 0.03436 307.59889 0.394 0.694
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.937
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular#boundary fit is singular
anova(lag_ins_age, lag_ins_age_randslope) #fit did not improve with inclusion of random slopes## Data: lag_dat
## Models:
## lag_ins_age: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_ins_age_randslope: Ins_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_ins_age 4 1013.3 1028.3 -502.64 1005.3
## lag_ins_age_randslope 6 1016.6 1039.2 -502.28 1004.6 0.7233 2 0.6965lag_glx_age<- lmer(Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_glx_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1396.1 1411.1 -694.0 1388.1 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6700 -0.6848 0.0127 0.6729 3.5338
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.000 0.000
## Residual 4.542 2.131
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 23.23676 0.39326 319.00000 59.088 < 2e-16 ***
## mri_age_y -0.25961 0.06309 319.00000 -4.115 4.93e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.953
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular#boundary fit is singular
lag_glx_age_randslope<- lmer(Glu_Gln_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_glx_age_randslope)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1399.5 1422.1 -693.7 1387.5 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6292 -0.6832 0.0209 0.6274 3.5685
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## subj_id (Intercept) 1.32593 1.1515
## mri_age_y 0.03495 0.1869 -1.00
## Residual 4.41427 2.1010
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 23.19469 0.41423 77.21723 55.994 < 2e-16 ***
## mri_age_y -0.24923 0.06668 67.41877 -3.738 0.000385 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.958
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular#boundary fit is singular
anova(lag_glx_age, lag_glx_age_randslope) #fit did not improve with inclusion of random slopes## Data: lag_dat
## Models:
## lag_glx_age: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_glx_age_randslope: Glu_Gln_conc_molal ~ mri_age_y + (1 + mri_age_y | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_glx_age 4 1396.1 1411.1 -694.04 1388.1
## lag_glx_age_randslope 6 1399.5 1422.1 -693.74 1387.5 0.5942 2 0.743##Tissue fraction analysis Calculate tissue fraction GM/(GM+WM) for each subject/voxel Add tissue fraction as fixed effect in lmer models ACG
library(ggplot2)
library(EnvStats)
library(lme4)
library(lmerTest)
#ACG
#read in data
acg_dat <- read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/acg_data_molal_outliers_removed_n113.csv")
#create tissue fraction variable
acg_dat$TissueFraction_GMWM <- acg_dat$fGM/(acg_dat$fGM + acg_dat$fWM)
#plot tissue fraction to check for issues
plot_ACG_tissuefraction <- ggplot(acg_dat, aes(TissueFraction_GMWM))
plot_ACG_tissuefraction+ geom_density() #big skew
plot_ACG_tissuefraction + geom_boxplot() #4 potential outliers
rosnerTest(acg_dat$TissueFraction_GMWM, k = 4)## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4
## 7.315531 5.204705 3.854675 3.544898
##
## $sample.size
## [1] 113
##
## $parameters
## k
## 4
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4
## 3.425263 3.422302 3.419309 3.416284
##
## $n.outliers
## [1] 4
##
## $alternative
## [1] "Up to 4 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 0.9390198 0.9661804 0.9655023 0.9366803 0.8741955 0.9311618 0.9320337
## [8] 0.9397518 0.9548737 0.9669866 0.9444564 0.9138117 0.9490174 0.9573432
## [15] 0.9497333 0.9748258 0.9941847 0.9486011 0.9602519 0.9612505 0.9732704
## [22] 0.9790679 0.9479635 0.8046696 0.9719504 0.9435072 0.9860912 0.9252084
## [29] 0.8946503 0.9275092 0.9413340 0.9619213 0.9344712 0.9996141 0.9611735
## [36] 0.9812598 0.6594113 0.9283294 0.9853063 0.9976770 0.9481254 0.9628117
## [43] 0.9347577 0.9444220 0.9509674 0.9248160 0.9468277 0.9750073 0.9871912
## [50] 0.9848082 0.9777584 0.9709937 0.9427788 0.9516720 0.9772194 0.8991193
## [57] 0.9724672 0.9533810 0.9848684 0.9562166 0.9898357 0.9339068 0.9541988
## [64] 0.9779450 0.9755787 0.9424898 0.9860622 0.9440049 0.9238023 0.9511135
## [71] 0.9519491 0.9895463 0.9500137 0.9847202 0.9509092 0.9630685 0.9736999
## [78] 0.9281849 0.9266655 0.9872865 0.9798646 0.9595847 0.9541502 0.9831991
## [85] 0.9377129 0.9309479 0.9847561 0.9729524 0.9754696 0.9392829 0.9671345
## [92] 0.9700358 0.8595117 0.9559970 0.9768332 0.9787093 0.9556623 0.9457820
## [99] 0.9936215 0.9881743 0.9804760 0.9764555 0.9591807 0.9290081 0.9073076
## [106] 0.9651388 0.9776632 0.9805040 0.9636278 0.9544220 0.9467918 0.9651846
## [113] 0.9882172
##
## $data.name
## [1] "acg_dat$TissueFraction_GMWM"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 0.9530867 0.04014410 0.6594113 37 7.315531 3.425263 TRUE
## 2 1 0.9557088 0.02901974 0.8046696 24 5.204705 3.422302 TRUE
## 3 2 0.9570695 0.02530896 0.8595117 93 3.854675 3.419309 TRUE
## 4 3 0.9579564 0.02362857 0.8741955 5 3.544898 3.416284 TRUE
##
## attr(,"class")
## [1] "gofOutlier"#4 outliers: Obs. Num 37 (PS14_090), 24 (PS14_050), 93 (PS18_1601), 5 (PS14_009) - double check cogregistration
#PS14_090 is shifted to the left and includes substantial WM - EXCLUDE
#PS14_050 is shifted to the left and T1 quality is poor - EXCLUDE
#PS18_1601 head tilted so voxel captures a bit of WM in the left hemi., a little anterior; doesn't seem terrible
#PS14_009 looks like chin was tipped down toward chest a bit, corner of voxel captures some CC white matter but seems ok
#save database with tissue fraction
write.csv(acg_dat, "/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/acg_data_molal_outliers_removed_n113_tissuefraction.csv", row.names = FALSE)
#exclude scans PS14_090 and PS14_050 due to outlier status coregistration issues
acg_dat_tf<-read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/acg_data_molal_outliers_removed_n113_tissuefraction.csv")
acg_dat_tf["37", "TissueFraction_GMWM"] <- NA
acg_dat_tf["24", "TissueFraction_GMWM"] <- NA
acg_dat_tf <- subset(acg_dat_tf, study_code != c("PS14_090", "PS14_050"))## Warning in study_code != c("PS14_090", "PS14_050"): longer object length is not
## a multiple of shorter object length#check correlation between tissue fraction and age
acg_TisFrac_age <- lmer(TissueFraction_GMWM ~ mri_age_y + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_TisFrac_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: TissueFraction_GMWM ~ mri_age_y + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## -499.0 -488.1 253.5 -507.0 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.98848 -0.32600 0.01794 0.42686 1.85387
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.0004829 0.02197
## Residual 0.0001820 0.01349
## Number of obs: 111, groups: subj_id, 100
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 9.552e-01 8.892e-03 1.002e+02 107.423 <2e-16 ***
## mri_age_y 3.986e-04 2.077e-03 9.572e+01 0.192 0.848
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.957#not sig.
#Add tissue fraction as fixed effect in ACG models
#NAA
#basic age model
acg_naa_age <- lmer(NAA_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_naa_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 364.2 375.0 -178.1 356.2 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.78812 -0.44480 0.04041 0.41080 1.88431
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.8872 0.9419
## Residual 0.6396 0.7998
## Number of obs: 111, groups: subj_id, 100
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 14.0978 0.4398 108.8460 32.052 <2e-16 ***
## mri_age_y 0.2044 0.1031 106.6834 1.982 0.0501 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.960#add tissue fraction
acg_naa_age_tf <- lmer(NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_naa_age_tf) #no sig. effect of tissue fraction## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 366.1 379.7 -178.1 356.1 106
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.78485 -0.44445 0.02803 0.41427 1.84724
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.9017 0.9496
## Residual 0.6272 0.7919
## Number of obs: 111, groups: subj_id, 100
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.1663 4.4412 110.4728 2.965 0.00371 **
## mri_age_y 0.2051 0.1030 106.1791 1.990 0.04915 *
## TissueFraction_GMWM 0.9703 4.6197 109.9872 0.210 0.83402
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.093
## TssFrc_GMWM -0.995 -0.002#compare models
anova(acg_naa_age, acg_naa_age_tf) #model fit not sig. different## Data: acg_dat_tf
## Models:
## acg_naa_age: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_naa_age_tf: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_naa_age 4 364.15 374.99 -178.08 356.15
## acg_naa_age_tf 5 366.11 379.66 -178.06 356.11 0.0424 1 0.8369#Cre
#basic age model
acg_cre_age <- lmer(Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_cre_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 324.5 335.4 -158.3 316.5 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.01085 -0.43475 0.01479 0.41135 1.81014
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.5452 0.7383
## Residual 0.5109 0.7147
## Number of obs: 111, groups: subj_id, 100
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 11.84959 0.36874 110.02636 32.136 <2e-16 ***
## mri_age_y 0.05989 0.08656 108.21331 0.692 0.49
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.961#add tissue fraction
acg_cre_age_tf <- lmer(Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_cre_age_tf)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 325.8 339.4 -157.9 315.8 106
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.05786 -0.42131 0.01284 0.42273 1.79300
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.5271 0.7260
## Residual 0.5203 0.7213
## Number of obs: 111, groups: subj_id, 100
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 14.88776 3.69640 109.55054 4.028 0.000104 ***
## mri_age_y 0.05874 0.08634 108.44592 0.680 0.497715
## TissueFraction_GMWM -3.17067 3.84260 108.86252 -0.825 0.411099
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.100
## TssFrc_GMWM -0.995 0.004#Cho
#basic age model
acg_cho_age <- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_cho_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 72.2 83.0 -32.1 64.2 106
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2025 -0.6306 -0.1413 0.6299 2.4918
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.02935 0.1713
## Residual 0.07690 0.2773
## Number of obs: 110, groups: subj_id, 99
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.76889 0.11894 109.88524 23.281 <2e-16 ***
## mri_age_y -0.06836 0.02799 109.96492 -2.442 0.0162 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.963#add tissue fraction
acg_cho_age_tf <- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_cho_age_tf)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 |
## subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 72.7 86.2 -31.4 62.7 105
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2564 -0.5970 -0.1487 0.6628 2.5884
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.02784 0.1668
## Residual 0.07690 0.2773
## Number of obs: 110, groups: subj_id, 99
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 4.19678 1.17742 109.25030 3.564 0.000542 ***
## mri_age_y -0.06876 0.02781 109.97149 -2.473 0.014935 *
## TissueFraction_GMWM -1.49055 1.22297 109.07759 -1.219 0.225551
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.109
## TssFrc_GMWM -0.995 0.013#compare models
anova(acg_cho_age, acg_cho_age_tf) #fit not significantly different## Data: acg_dat_tf
## Models:
## acg_cho_age: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## acg_cho_age_tf: Cho_GPC_PCh_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## acg_cho_age 4 72.191 82.993 -32.096 64.191
## acg_cho_age_tf 5 72.718 86.221 -31.359 62.718 1.473 1 0.2249#Ins
#basic age model
acg_ins_age <- lmer(Ins_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_ins_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 407.8 418.6 -199.9 399.8 106
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5246 -0.4721 0.0711 0.5130 1.5756
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.9237 0.9611
## Residual 1.3525 1.1630
## Number of obs: 110, groups: subj_id, 99
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 8.5924 0.5484 109.9754 15.667 <2e-16 ***
## mri_age_y -0.2258 0.1294 109.4869 -1.746 0.0837 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.962#add tissue fraction
acg_ins_age_tf <- lmer(Ins_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_ins_age_tf)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 409.4 422.9 -199.7 399.4 105
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.52269 -0.47942 0.08856 0.57698 1.49140
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.9414 0.9702
## Residual 1.3276 1.1522
## Number of obs: 110, groups: subj_id, 99
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 4.8246 5.4937 108.8944 0.878 0.3818
## mri_age_y -0.2249 0.1291 109.3760 -1.742 0.0844 .
## TissueFraction_GMWM 3.9350 5.7034 108.5006 0.690 0.4917
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.116
## TssFrc_GMWM -0.995 0.021#Glx
#basic age model
acg_glx_age <- lmer(Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_glx_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 530.3 541.1 -261.1 522.3 107
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.48912 -0.37848 -0.00398 0.44809 1.68844
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 4.601 2.145
## Residual 2.344 1.531
## Number of obs: 111, groups: subj_id, 100
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 25.9373 0.9245 105.2183 28.055 <2e-16 ***
## mri_age_y 0.2288 0.2164 101.8687 1.057 0.293
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.959#add tissue fraction
acg_glx_age_tf <- lmer(Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=acg_dat_tf, REML=FALSE)
summary(acg_glx_age_tf)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: acg_dat_tf
##
## AIC BIC logLik deviance df.resid
## 530.5 544.0 -260.2 520.5 106
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.43181 -0.34955 -0.03061 0.41546 1.49601
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 4.837 2.199
## Residual 2.068 1.438
## Number of obs: 111, groups: subj_id, 100
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 13.2424 9.3236 110.9911 1.420 0.158
## mri_age_y 0.2289 0.2135 98.4861 1.072 0.286
## TissueFraction_GMWM 13.2676 9.7105 110.9341 1.366 0.175
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.080
## TssFrc_GMWM -0.995 -0.014LAG Tissue Fraction models
#LAG
#read in data
lag_dat<-read.csv("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/lag_data_molal_outliers_removed_n320.csv")
#create tissue fraction variable
lag_dat$TissueFraction_GMWM <- lag_dat$fGM/(lag_dat$fGM + lag_dat$fWM)
#plot tissue fraction to check for issues
plot_LAG_tissuefraction <- ggplot(lag_dat, aes(TissueFraction_GMWM))
plot_LAG_tissuefraction+ geom_density() #normal
plot_LAG_tissuefraction + geom_boxplot() #6 potential outliers
rosnerTest(lag_dat$TissueFraction_GMWM, k = 6) #no significant outliers## $distribution
## [1] "Normal"
##
## $statistic
## R.1 R.2 R.3 R.4 R.5 R.6
## 3.280842 3.075086 2.985173 2.935103 2.963794 2.932421
##
## $sample.size
## [1] 321
##
## $parameters
## k
## 6
##
## $alpha
## [1] 0.05
##
## $crit.value
## lambda.1 lambda.2 lambda.3 lambda.4 lambda.5 lambda.6
## 3.742579 3.741706 3.740829 3.739949 3.739067 3.738181
##
## $n.outliers
## [1] 0
##
## $alternative
## [1] "Up to 6 observations are not\n from the same Distribution."
##
## $method
## [1] "Rosner's Test for Outliers"
##
## $data
## [1] 0.7657305 0.6515751 0.5744539 0.6894155 0.8381194 0.5951315 0.7008113
## [8] 0.6813746 0.5352570 0.6558286 0.7097519 0.7819635 0.7021001 0.7039674
## [15] 0.4764324 0.6028844 0.6327450 0.6087942 0.6968812 0.8026664 0.7009372
## [22] 0.7208761 0.6605230 0.7324027 0.5453879 0.7620629 0.6701616 0.6375845
## [29] 0.5632865 0.6708753 0.6927008 0.7798970 0.7618192 0.6477681 0.6613518
## [36] 0.5340622 0.5942859 0.6943742 0.7009750 0.5879671 0.4411657 0.6156132
## [43] 0.6142354 0.7259517 0.5979064 0.5685179 0.5271897 0.5355836 0.7237169
## [50] 0.7571831 0.5285119 0.6751229 0.5119628 0.7521085 0.6615893 0.5425676
## [57] 0.6627217 0.6549878 0.7279301 0.7256598 0.7477565 0.6395054 0.7213695
## [64] 0.7384663 0.7633997 0.7529087 0.8002460 0.8145222 0.7366533 0.7977780
## [71] 0.6151766 0.7154834 0.7754343 0.8142228 0.7207045 0.6896770 0.5855033
## [78] 0.4519376 0.3757102 0.6561017 0.7878613 0.7377081 0.8047878 0.7816798
## [85] 0.8279033 0.6753913 0.7666730 0.5795116 0.7313181 0.6599876 0.7201454
## [92] 0.6696166 0.5834821 0.6485976 0.5740400 0.7069867 0.7142496 0.6394643
## [99] 0.7231628 0.6282001 0.7420446 0.5742007 0.7048192 0.8067634 0.8402937
## [106] 0.6709200 0.7285484 0.6456249 0.7123905 0.7182381 0.6336725 0.7074584
## [113] 0.5114034 0.7548622 0.7512744 0.6852121 0.6741392 0.4619200 0.5944252
## [120] 0.6890812 0.5371624 0.6516809 0.6145975 0.6753211 0.5748713 0.7172832
## [127] 0.7275696 0.4961396 0.6162007 0.4003669 0.6812224 0.6283159 0.6387185
## [134] 0.7522739 0.4447576 0.5459764 0.5522673 0.5951632 0.5682089 0.5371529
## [141] 0.6615547 0.6791234 0.7036370 0.6886679 0.6631862 0.6341630 0.6738518
## [148] 0.7549798 0.6955773 0.6542445 0.6925320 0.6505739 0.6425831 0.6677008
## [155] 0.7839529 0.6553603 0.5432005 0.6122120 0.6548850 0.6367154 0.7679106
## [162] 0.7485923 0.7184663 0.6312047 0.7349290 0.7043475 0.4988189 0.8401854
## [169] 0.6577281 0.8438765 0.8390803 0.6282136 0.6405363 0.7191654 0.6334385
## [176] 0.7081752 0.6659769 0.6107064 0.6347357 0.7290806 0.7344315 0.6660147
## [183] 0.6251258 0.5997631 0.5785047 0.7984798 0.5070360 0.6098369 0.6299690
## [190] 0.6516707 0.5494253 0.6308773 0.7664214 0.5915382 0.5878598 0.8187175
## [197] 0.6741752 0.7564550 0.6443908 0.7395504 0.5204329 0.6170451 0.7002564
## [204] 0.6246470 0.6681476 0.5771348 0.7295070 0.7423764 0.6384295 0.5430082
## [211] 0.6707773 0.6881614 0.7608897 0.5543782 0.6804216 0.5896538 0.6427618
## [218] 0.6536086 0.7261262 0.6317773 0.6652328 0.6605223 0.5295042 0.6174694
## [225] 0.6425888 0.5149187 0.7329097 0.5653457 0.5104815 0.6635529 0.7475020
## [232] 0.7604907 0.6907149 0.5576525 0.6882708 0.7095940 0.5975811 0.6217704
## [239] 0.6017863 0.6938276 0.5893048 0.4639501 0.7554426 0.5394351 0.5730658
## [246] 0.5198421 0.6217217 0.5560520 0.7715035 0.7450276 0.7578195 0.7198068
## [253] 0.6385403 0.4769133 0.6216745 0.7215428 0.7135267 0.5015701 0.5463705
## [260] 0.3420581 0.4840778 0.4052625 0.3302661 0.3499400 0.6039142 0.4129804
## [267] 0.5315517 0.5809198 0.4792521 0.5046752 0.5684662 0.5218814 0.3154764
## [274] 0.3403429 0.5024469 0.3751392 0.6127580 0.7315085 0.4137321 0.4513800
## [281] 0.5718875 0.7465811 0.5551458 0.4922670 0.5978026 0.6036407 0.7885750
## [288] 0.7143374 0.5734673 0.7879821 0.5626092 0.6622175 0.4524319 0.4110109
## [295] 0.6597638 0.5684290 0.4234816 0.4456416 0.5029994 0.5793653 0.6290320
## [302] 0.4115474 0.6699743 0.7957488 0.7141763 0.7477478 0.5862060 0.7635972
## [309] 0.7800932 0.5786560 0.5839172 0.6023700 0.5311156 0.5390859 0.6565085
## [316] 0.8540765 0.7585903 0.5580707 0.6451797 0.4621925 0.9951620
##
## $data.name
## [1] "lag_dat$TissueFraction_GMWM"
##
## $bad.obs
## [1] 0
##
## $all.stats
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 0.6422445 0.1075692 0.9951620 321 3.280842 3.742579 FALSE
## 2 1 0.6411416 0.1059044 0.3154764 273 3.075086 3.741706 FALSE
## 3 2 0.6421625 0.1044818 0.3302661 263 2.985173 3.740829 FALSE
## 4 3 0.6431433 0.1031652 0.3403429 274 2.935103 3.739949 FALSE
## 5 4 0.6440985 0.1019101 0.3420581 260 2.963794 3.739067 FALSE
## 6 5 0.6450543 0.1006385 0.3499400 264 2.932421 3.738181 FALSE
##
## attr(,"class")
## [1] "gofOutlier"write.csv(lag_dat, "/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/data/lag_data_molal_outliers_removed_n320_tissuefraction.csv")
#check correlation between tissue fraction and age
lag_TisFrac_age <- lmer(TissueFraction_GMWM ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_TisFrac_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: TissueFraction_GMWM ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## -562.8 -547.7 285.4 -570.8 317
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.7967 -0.6352 0.0458 0.6787 3.5935
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.0006557 0.02561
## Residual 0.0093041 0.09646
## Number of obs: 321, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.770798 0.018475 306.068797 41.722 < 2e-16 ***
## mri_age_y -0.021705 0.002927 320.934582 -7.417 1.08e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.943#tissue fraction is significantly negatively associated with age
#in a simple correlation test
cor.test(lag_dat$mri_age_y, lag_dat$TissueFraction_GMWM) # also significant: r = -.37##
## Pearson's product-moment correlation
##
## data: lag_dat$mri_age_y and lag_dat$TissueFraction_GMWM
## t = -7.1217, df = 319, p-value = 7.113e-12
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.4611519 -0.2719374
## sample estimates:
## cor
## -0.3703805#Add tissue fraction as a fixed effect in LAG models
#basic age model NAA
lag_naa_age<- lmer(NAA_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
#add tissue fraction
lag_naa_age_TisFrac <- lmer(NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_naa_age_TisFrac)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 858.0 876.8 -424.0 848.0 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9925 -0.5718 -0.0199 0.5594 3.9440
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1630 0.4037
## Residual 0.7184 0.8476
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 11.22872 0.43459 317.97211 25.838 < 2e-16 ***
## mri_age_y 0.17233 0.02961 317.97750 5.819 1.44e-08 ***
## TissueFraction_GMWM 3.35162 0.51255 313.38178 6.539 2.51e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.703
## TssFrc_GMWM -0.915 0.394#compare model fit
anova(lag_naa_age, lag_naa_age_TisFrac) #fit significantly improved## Data: lag_dat
## Models:
## lag_naa_age: NAA_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_naa_age_TisFrac: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_naa_age 4 895.96 911.01 -443.98 887.96
## lag_naa_age_TisFrac 5 857.96 876.77 -423.98 847.96 39.996 1 2.545e-10
##
## lag_naa_age
## lag_naa_age_TisFrac ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#add tissue fraction interaction with age
lag_naa_age_TisFracXage <- lmer(NAA_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_naa_age_TisFracXage)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 857.3 879.9 -422.6 845.3 312
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1567 -0.5361 -0.0447 0.5566 4.0061
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1599 0.3999
## Residual 0.7134 0.8447
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 9.7907 0.9790 313.2768 10.001 < 2e-16 ***
## mri_age_y 0.3877 0.1348 307.1007 2.877 0.004302 **
## TissueFraction_GMWM 5.7025 1.5239 311.3849 3.742 0.000217 ***
## mri_age_y:TissueFraction_GMWM -0.3606 0.2202 305.3472 -1.638 0.102532
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.943
## TssFrc_GMWM -0.981 0.948
## m__:TF_GMWM 0.897 -0.976 -0.942#interaction term for agexTissueFraction not sig.
#compare model fit
anova(lag_naa_age_TisFrac, lag_naa_age_TisFracXage) #not sig., stick with age + TissueFraction model## Data: lag_dat
## Models:
## lag_naa_age_TisFrac: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## lag_naa_age_TisFracXage: NAA_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_naa_age_TisFrac 5 857.96 876.77 -423.98 847.96
## lag_naa_age_TisFracXage 6 857.29 879.86 -422.65 845.29 2.6694 1
## Pr(>Chisq)
## lag_naa_age_TisFrac
## lag_naa_age_TisFracXage 0.1023#add sex to tissue fraction model
lag_naa_age_TisFrac_sex <- lmer(NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_naa_age_TisFrac_sex) #sex not significant## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + female + (1 |
## subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 859.9 882.4 -423.9 847.9 312
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9834 -0.5775 -0.0126 0.5684 3.9313
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1616 0.4020
## Residual 0.7190 0.8479
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 11.24098 0.43642 317.96672 25.757 < 2e-16 ***
## mri_age_y 0.17281 0.02965 317.99247 5.828 1.37e-08 ***
## TissueFraction_GMWM 3.35735 0.51294 313.49233 6.545 2.42e-10 ***
## female -0.03974 0.13192 97.72196 -0.301 0.764
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.695
## TssFrc_GMWM -0.907 0.395
## female -0.093 -0.052 -0.039#compare model fit
anova(lag_naa_age_TisFrac, lag_naa_age_TisFrac_sex) #not sig., stick with age + TissueFraction model## Data: lag_dat
## Models:
## lag_naa_age_TisFrac: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## lag_naa_age_TisFrac_sex: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_naa_age_TisFrac 5 857.96 876.77 -423.98 847.96
## lag_naa_age_TisFrac_sex 6 859.87 882.44 -423.94 847.87 0.0901 1
## Pr(>Chisq)
## lag_naa_age_TisFrac
## lag_naa_age_TisFrac_sex 0.7641#Cre
#basic age model
lag_cre_age<- lmer(Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cre_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 848.2 863.3 -420.1 840.2 316
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9764 -0.5307 -0.1261 0.5525 4.1548
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.08477 0.2912
## Residual 0.73683 0.8584
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 10.49912 0.16805 308.49608 62.476 <2e-16 ***
## mri_age_y -0.04183 0.02646 319.87112 -1.581 0.115
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.937#add tissue fraction
lag_cre_age_TisFrac<- lmer(Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cre_age_TisFrac)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 760.2 779.0 -375.1 750.2 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7493 -0.5056 0.0270 0.4419 4.8839
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.08547 0.2923
## Residual 0.54120 0.7357
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 7.04677 0.36923 319.70413 19.085 <2e-16 ***
## mri_age_y 0.05661 0.02497 319.62390 2.267 0.0241 *
## TissueFraction_GMWM 4.47115 0.43745 317.80861 10.221 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.701
## TssFrc_GMWM -0.918 0.393#compare model fit
anova(lag_cre_age, lag_cre_age_TisFrac) #sig.## Data: lag_dat
## Models:
## lag_cre_age: Cre_PCr_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_cre_age_TisFrac: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_cre_age 4 848.23 863.31 -420.12 840.23
## lag_cre_age_TisFrac 5 760.20 779.04 -375.10 750.20 90.029 1 < 2.2e-16
##
## lag_cre_age
## lag_cre_age_TisFrac ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#add interaction
lag_cre_age_TisFrac_int<- lmer(Cre_PCr_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cre_age_TisFrac_int)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 761.9 784.5 -375.0 749.9 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7360 -0.5018 0.0276 0.4507 4.9009
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.08463 0.2909
## Residual 0.54125 0.7357
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.65386 0.83511 317.57482 7.968 2.93e-14
## mri_age_y 0.11569 0.11532 312.49738 1.003 0.316524
## TissueFraction_GMWM 5.11282 1.29936 316.10788 3.935 0.000102
## mri_age_y:TissueFraction_GMWM -0.09884 0.18836 310.83144 -0.525 0.600131
##
## (Intercept) ***
## mri_age_y
## TissueFraction_GMWM ***
## mri_age_y:TissueFraction_GMWM
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.943
## TssFrc_GMWM -0.981 0.948
## m__:TF_GMWM 0.897 -0.976 -0.942#compare model fit
anova(lag_cre_age_TisFrac, lag_cre_age_TisFrac_int) #not sig.## Data: lag_dat
## Models:
## lag_cre_age_TisFrac: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## lag_cre_age_TisFrac_int: Cre_PCr_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_cre_age_TisFrac 5 760.20 779.04 -375.10 750.20
## lag_cre_age_TisFrac_int 6 761.93 784.54 -374.96 749.93 0.2748 1
## Pr(>Chisq)
## lag_cre_age_TisFrac
## lag_cre_age_TisFrac_int 0.6002#add sex to tissue fraction model
lag_cre_age_TisFrac_sex<- lmer(Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cre_age_TisFrac_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + female +
## (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 762.0 784.6 -375.0 750.0 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7681 -0.5123 0.0220 0.4570 4.9061
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.0852 0.2919
## Residual 0.5410 0.7355
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 7.03303 0.37047 319.29640 18.984 <2e-16 ***
## mri_age_y 0.05601 0.02501 319.56737 2.240 0.0258 *
## TissueFraction_GMWM 4.46392 0.43763 317.89203 10.200 <2e-16 ***
## female 0.04644 0.10629 90.75521 0.437 0.6632
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.693
## TssFrc_GMWM -0.911 0.394
## female -0.085 -0.056 -0.037#compare model fit
anova(lag_cre_age_TisFrac, lag_cre_age_TisFrac_sex) #not sig.## Data: lag_dat
## Models:
## lag_cre_age_TisFrac: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## lag_cre_age_TisFrac_sex: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_cre_age_TisFrac 5 760.20 779.04 -375.10 750.20
## lag_cre_age_TisFrac_sex 6 762.01 784.62 -375.01 750.01 0.1908 1
## Pr(>Chisq)
## lag_cre_age_TisFrac
## lag_cre_age_TisFrac_sex 0.6622#Cho
#basic age model
lag_cho_age<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 51.3 66.4 -21.7 43.3 316
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6539 -0.5810 -0.0562 0.4659 4.1202
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.0149 0.1221
## Residual 0.0561 0.2368
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.346670 0.049038 310.144986 47.854 < 2e-16 ***
## mri_age_y -0.023315 0.007581 315.773363 -3.075 0.00229 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.921#add tissue fraction
lag_cho_age_TisFrac<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_age_TisFrac)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 |
## subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 25.3 44.2 -7.7 15.3 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5688 -0.6208 -0.0348 0.5176 4.6606
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.01202 0.1096
## Residual 0.05232 0.2287
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.929074 0.117156 319.933789 25.002 < 2e-16 ***
## mri_age_y -0.040403 0.007925 319.967590 -5.098 5.89e-07 ***
## TissueFraction_GMWM -0.749693 0.138214 313.218002 -5.424 1.17e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.704
## TssFrc_GMWM -0.917 0.399#compare model fit
anova(lag_cho_age, lag_cho_age_TisFrac) #sig## Data: lag_dat
## Models:
## lag_cho_age: Cho_GPC_PCh_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_cho_age_TisFrac: Cho_GPC_PCh_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_cho_age 4 51.345 66.419 -21.6727 43.345
## lag_cho_age_TisFrac 5 25.333 44.174 -7.6663 15.333 28.013 1 1.205e-07
##
## lag_cho_age
## lag_cho_age_TisFrac ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#add interaction
lag_cho_age_TisFrac_int<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_age_TisFrac_int) #interaction is sig. (p=.046)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 |
## subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 23.3 46.0 -5.7 11.3 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5480 -0.6259 -0.0670 0.5067 4.4339
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.01168 0.1081
## Residual 0.05179 0.2276
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.45851 0.26231 312.15748 9.373 <2e-16
## mri_age_y 0.03027 0.03614 303.65993 0.838 0.403
## TissueFraction_GMWM 0.01894 0.40778 309.31482 0.046 0.963
## mri_age_y:TissueFraction_GMWM -0.11820 0.05899 301.29296 -2.004 0.046
##
## (Intercept) ***
## mri_age_y
## TissueFraction_GMWM
## mri_age_y:TissueFraction_GMWM *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.943
## TssFrc_GMWM -0.981 0.948
## m__:TF_GMWM 0.896 -0.976 -0.942#compare model fit
anova(lag_cho_age_TisFrac, lag_cho_age_TisFrac_int) #sig. (p<.046)## Data: lag_dat
## Models:
## lag_cho_age_TisFrac: Cho_GPC_PCh_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## lag_cho_age_TisFrac_int: Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_cho_age_TisFrac 5 25.333 44.174 -7.6663 15.333
## lag_cho_age_TisFrac_int 6 23.345 45.955 -5.6724 11.345 3.9878 1
## Pr(>Chisq)
## lag_cho_age_TisFrac
## lag_cho_age_TisFrac_int 0.04583 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#add sex to model
lag_cho_age_TisFrac_int_sex<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_cho_age_TisFrac_int_sex) #not sig.## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + female +
## (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 24.2 50.6 -5.1 10.2 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5865 -0.6287 -0.0611 0.5043 4.4825
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.01157 0.1076
## Residual 0.05164 0.2273
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.485820 0.263137 312.684444 9.447 <2e-16
## mri_age_y 0.028390 0.036124 303.303507 0.786 0.4325
## TissueFraction_GMWM -0.001475 0.407568 309.052784 -0.004 0.9971
## female -0.037481 0.035376 77.480663 -1.060 0.2927
## mri_age_y:TissueFraction_GMWM -0.114251 0.059018 300.819739 -1.936 0.0538
##
## (Intercept) ***
## mri_age_y
## TissueFraction_GMWM
## female
## mri_age_y:TissueFraction_GMWM .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW female
## mri_age_y -0.942
## TssFrc_GMWM -0.980 0.948
## female -0.098 0.050 0.047
## m__:TF_GMWM 0.896 -0.976 -0.942 -0.064#compare model fit
anova(lag_cho_age_TisFrac_int, lag_cho_age_TisFrac_int_sex) #not sig.## Data: lag_dat
## Models:
## lag_cho_age_TisFrac_int: Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## lag_cho_age_TisFrac_int_sex: Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_cho_age_TisFrac_int 6 23.345 45.955 -5.6724 11.345
## lag_cho_age_TisFrac_int_sex 7 24.224 50.603 -5.1122 10.225 1.1204 1
## Pr(>Chisq)
## lag_cho_age_TisFrac_int
## lag_cho_age_TisFrac_int_sex 0.2898#Ins
#basic age model
lag_ins_age<- lmer(Ins_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_ins_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1013.3 1028.3 -502.6 1005.3 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3984 -0.5594 -0.0111 0.5911 4.6516
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1768 0.4205
## Residual 1.2227 1.1057
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.28473 0.21938 307.28604 28.648 <2e-16 ***
## mri_age_y 0.01632 0.03444 318.45543 0.474 0.636
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.934#add tissue fraction
lag_ins_age_TisFrac<- lmer(Ins_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_ins_age_TisFrac) #tissue fraction sig.; age not sig.## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1005.5 1024.3 -497.7 995.5 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3841 -0.5608 -0.0382 0.5912 4.7226
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1568 0.396
## Residual 1.1958 1.094
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 4.70590 0.54440 318.16597 8.644 2.69e-16 ***
## mri_age_y 0.06153 0.03685 318.01860 1.670 0.09596 .
## TissueFraction_GMWM 2.04247 0.64587 317.48996 3.162 0.00172 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.700
## TssFrc_GMWM -0.918 0.391#compare model fit
anova(lag_ins_age, lag_ins_age_TisFrac) #sig.## Data: lag_dat
## Models:
## lag_ins_age: Ins_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_ins_age_TisFrac: Ins_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_ins_age 4 1013.3 1028.3 -502.64 1005.28
## lag_ins_age_TisFrac 5 1005.5 1024.3 -497.73 995.47 9.812 1 0.001734 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#add interaction
lag_ins_age_TisFrac_int<- lmer(Ins_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_ins_age_TisFrac_int) #tissue fraction, age & interaction all sig.## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1001.6 1024.2 -494.8 989.6 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3911 -0.5979 -0.0053 0.5734 4.7179
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1767 0.4203
## Residual 1.1584 1.0763
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.0286 1.2206 316.2405 1.662 0.097514 .
## mri_age_y 0.4654 0.1686 310.2011 2.761 0.006114 **
## TissueFraction_GMWM 6.4126 1.8995 314.4820 3.376 0.000828 ***
## mri_age_y:TissueFraction_GMWM -0.6751 0.2755 308.0776 -2.450 0.014822 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.943
## TssFrc_GMWM -0.981 0.948
## m__:TF_GMWM 0.897 -0.976 -0.942#compare model fit
anova(lag_ins_age_TisFrac, lag_ins_age_TisFrac_int) #sig.## Data: lag_dat
## Models:
## lag_ins_age_TisFrac: Ins_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## lag_ins_age_TisFrac_int: Ins_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_ins_age_TisFrac 5 1005.5 1024.3 -497.73 995.47
## lag_ins_age_TisFrac_int 6 1001.6 1024.2 -494.80 989.60 5.8652 1
## Pr(>Chisq)
## lag_ins_age_TisFrac
## lag_ins_age_TisFrac_int 0.01544 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#add sex
lag_ins_age_TisFrac_int_sex<- lmer(Ins_conc_molal ~ mri_age_y * TissueFraction_GMWM + female + (1 | subj_id), data=lag_dat, REML=FALSE)
summary(lag_ins_age_TisFrac_int_sex) #tissue fraction, age & interaction all sig.## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Ins_conc_molal ~ mri_age_y * TissueFraction_GMWM + female + (1 |
## subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1003.6 1030.0 -494.8 989.6 312
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3882 -0.5982 -0.0079 0.5748 4.7194
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1772 0.4209
## Residual 1.1581 1.0761
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.023474 1.226394 316.309945 1.650 0.099948
## mri_age_y 0.465876 0.168831 309.440740 2.759 0.006135
## TissueFraction_GMWM 6.417150 1.901892 314.005179 3.374 0.000833
## female 0.005839 0.155173 78.119816 0.038 0.970081
## mri_age_y:TissueFraction_GMWM -0.675931 0.276096 307.149063 -2.448 0.014918
##
## (Intercept) .
## mri_age_y **
## TissueFraction_GMWM ***
## female
## mri_age_y:TissueFraction_GMWM *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW female
## mri_age_y -0.942
## TssFrc_GMWM -0.980 0.948
## female -0.097 0.054 0.051
## m__:TF_GMWM 0.897 -0.976 -0.942 -0.068#compare model fit
anova(lag_ins_age_TisFrac_int, lag_ins_age_TisFrac_int_sex) #not sig.## Data: lag_dat
## Models:
## lag_ins_age_TisFrac_int: Ins_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## lag_ins_age_TisFrac_int_sex: Ins_conc_molal ~ mri_age_y * TissueFraction_GMWM + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_ins_age_TisFrac_int 6 1001.6 1024.2 -494.8 989.6
## lag_ins_age_TisFrac_int_sex 7 1003.6 1030.0 -494.8 989.6 0.0014 1
## Pr(>Chisq)
## lag_ins_age_TisFrac_int
## lag_ins_age_TisFrac_int_sex 0.9704#Glx
#basic age model
lag_glx_age<- lmer(Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_glx_age)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1396.1 1411.1 -694.0 1388.1 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6700 -0.6848 0.0127 0.6729 3.5338
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.000 0.000
## Residual 4.542 2.131
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 23.23676 0.39326 319.00000 59.088 < 2e-16 ***
## mri_age_y -0.25961 0.06309 319.00000 -4.115 4.93e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## mri_age_y -0.953
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular#add tissue fraction
lag_glx_age_TisFrac<- lmer(Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_glx_age_TisFrac)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1238.1 1257.0 -614.1 1228.1 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5542 -0.6325 -0.0043 0.6145 4.7331
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.000 0.000
## Residual 2.751 1.659
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 12.97325 0.77521 319.00000 16.735 <2e-16 ***
## mri_age_y 0.02096 0.05282 319.00000 0.397 0.692
## TissueFraction_GMWM 13.37772 0.92834 319.00000 14.410 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.688
## TssFrc_GMWM -0.919 0.369
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular#compare models
anova(lag_glx_age, lag_glx_age_TisFrac) #sig## Data: lag_dat
## Models:
## lag_glx_age: Glu_Gln_conc_molal ~ mri_age_y + (1 | subj_id)
## lag_glx_age_TisFrac: Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## lag_glx_age 4 1396.1 1411.1 -694.04 1388.1
## lag_glx_age_TisFrac 5 1238.1 1257.0 -614.07 1228.1 159.94 1 < 2.2e-16
##
## lag_glx_age
## lag_glx_age_TisFrac ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#add interaction
lag_glx_age_TisFrac_int<- lmer(Glu_Gln_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_glx_age_TisFrac_int)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1237.7 1260.3 -612.8 1225.7 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7571 -0.6187 0.0099 0.5924 4.7144
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.00 0.000
## Residual 2.73 1.652
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 15.4568 1.7652 319.0000 8.756 < 2e-16 ***
## mri_age_y -0.3547 0.2458 319.0000 -1.443 0.149993
## TissueFraction_GMWM 9.3213 2.7527 319.0000 3.386 0.000797 ***
## mri_age_y:TissueFraction_GMWM 0.6292 0.4022 319.0000 1.565 0.118677
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.943
## TssFrc_GMWM -0.982 0.947
## m__:TF_GMWM 0.899 -0.977 -0.942
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular#compare models
anova(lag_glx_age_TisFrac, lag_glx_age_TisFrac_int) #not sig## Data: lag_dat
## Models:
## lag_glx_age_TisFrac: Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## lag_glx_age_TisFrac_int: Glu_Gln_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_glx_age_TisFrac 5 1238.1 1257.0 -614.07 1228.1
## lag_glx_age_TisFrac_int 6 1237.7 1260.3 -612.85 1225.7 2.4385 1
## Pr(>Chisq)
## lag_glx_age_TisFrac
## lag_glx_age_TisFrac_int 0.1184#add sex
lag_glx_age_TisFrac_sex<- lmer(Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + female + (1 | subj_id), data=lag_dat, REML=FALSE)## boundary (singular) fit: see ?isSingularsummary(lag_glx_age_TisFrac_sex)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + female +
## (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 1239.1 1261.7 -613.5 1227.1 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.4999 -0.6414 0.0131 0.6445 4.7946
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.000 0.000
## Residual 2.742 1.656
## Number of obs: 319, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 12.91764 0.77583 319.00000 16.650 <2e-16 ***
## mri_age_y 0.01866 0.05278 319.00000 0.353 0.724
## TissueFraction_GMWM 13.34681 0.92730 319.00000 14.393 <2e-16 ***
## female 0.19087 0.18604 319.00000 1.026 0.306
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.683
## TssFrc_GMWM -0.914 0.369
## female -0.070 -0.042 -0.032
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular#compare models
anova(lag_glx_age_TisFrac, lag_glx_age_TisFrac_sex) #not sig## Data: lag_dat
## Models:
## lag_glx_age_TisFrac: Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## lag_glx_age_TisFrac_sex: Glu_Gln_conc_molal ~ mri_age_y + TissueFraction_GMWM + female + (1 | subj_id)
## npar AIC BIC logLik deviance Chisq Df
## lag_glx_age_TisFrac 5 1238.1 1257.0 -614.07 1228.1
## lag_glx_age_TisFrac_sex 6 1239.1 1261.7 -613.54 1227.1 1.0509 1
## Pr(>Chisq)
## lag_glx_age_TisFrac
## lag_glx_age_TisFrac_sex 0.3053#check models for multicollinearity
library(performance)
check_collinearity(lag_naa_age_TisFrac) #Low## # Check for Multicollinearity
##
## Low Correlation
##
## Term VIF Increased SE Tolerance
## mri_age_y 1.18 1.09 0.84
## TissueFraction_GMWM 1.18 1.09 0.84check_collinearity(lag_cre_age_TisFrac) #Low## # Check for Multicollinearity
##
## Low Correlation
##
## Term VIF Increased SE Tolerance
## mri_age_y 1.18 1.09 0.85
## TissueFraction_GMWM 1.18 1.09 0.85check_collinearity(lag_cho_age_TisFrac_int) #interaction term causes inflated VIFs - check model without interaction## Warning: Model has interaction terms. VIFs might be inflated. You may check
## multicollinearity among predictors of a model without interaction terms.## # Check for Multicollinearity
##
## High Correlation
##
## Term VIF Increased SE Tolerance
## mri_age_y 25.02 5.00 0.04
## TissueFraction_GMWM 10.47 3.24 0.10
## mri_age_y:TissueFraction_GMWM 22.26 4.72 0.04check_collinearity(lag_cho_age_TisFrac) #Low## # Check for Multicollinearity
##
## Low Correlation
##
## Term VIF Increased SE Tolerance
## mri_age_y 1.19 1.09 0.84
## TissueFraction_GMWM 1.19 1.09 0.84check_collinearity(lag_ins_age_TisFrac) #Low## # Check for Multicollinearity
##
## Low Correlation
##
## Term VIF Increased SE Tolerance
## mri_age_y 1.18 1.09 0.85
## TissueFraction_GMWM 1.18 1.09 0.85check_collinearity(lag_glx_age_TisFrac) #Low## # Check for Multicollinearity
##
## Low Correlation
##
## Term VIF Increased SE Tolerance
## mri_age_y 1.16 1.08 0.86
## TissueFraction_GMWM 1.16 1.08 0.86##Unpack interactions Use sim_slopes function from interactions package to conduct simple slopes analysis on continuous x continuous interactions
library(interactions)
#LAG Cho: age x tissue fraction interaction
lag_cho_age_TisFrac_int<- lmer(Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
#run simple slopes analysis with Johnson-Neyman plot
sim_slopes(lag_cho_age_TisFrac_int, pred = mri_age_y, modx = TissueFraction_GMWM, jnplot = TRUE)## JOHNSON-NEYMAN INTERVAL
##
## When TissueFraction_GMWM is OUTSIDE the interval [-18.37, 0.45], the slope
## of mri_age_y is p < .05.
##
## Note: The range of observed values of TissueFraction_GMWM is [0.32, 1.00]
## SIMPLE SLOPES ANALYSIS
##
## Slope of mri_age_y when TissueFraction_GMWM = 0.5347131 (- 1 SD):
##
## Est. S.E. t val. p
## ------- ------ -------- ------
## -0.03 0.01 -3.78 0.00
##
## Slope of mri_age_y when TissueFraction_GMWM = 0.6424099 (Mean):
##
## Est. S.E. t val. p
## ------- ------ -------- ------
## -0.05 0.01 -5.50 0.00
##
## Slope of mri_age_y when TissueFraction_GMWM = 0.7501066 (+ 1 SD):
##
## Est. S.E. t val. p
## ------- ------ -------- ------
## -0.06 0.01 -4.89 0.00#save Johnson-Neyman plot
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/lag_cho_age_tf_int_JNplot.png")## Saving 7 x 5 in image#create interaction plot (make it pretty later)
#basic plot
Cho_age_tf_int_plot <- interact_plot(lag_cho_age_TisFrac_int, pred = mri_age_y, modx = TissueFraction_GMWM)
Cho_age_tf_int_plot
#linearity check
Cho_age_tf_int_plot_lin <- interact_plot(lag_cho_age_TisFrac_int, pred = mri_age_y, modx = TissueFraction_GMWM, modx.values="terciles", plot.points=TRUE, linearity.check=TRUE)## Medians of each tercile of TissueFraction_GMWM are 0.539, 0.655, 0.747Cho_age_tf_int_plot_lin## `geom_smooth()` using formula 'y ~ x'
#LAG Ins: age x tissue fraction interaction
lag_ins_age_TisFrac_int<- lmer(Ins_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 | subj_id), data=lag_dat, REML=FALSE)
#run simple slopes analysis with Johnson-Neyman plot
sim_slopes(lag_ins_age_TisFrac_int, pred = mri_age_y, modx = TissueFraction_GMWM, jnplot = TRUE)## JOHNSON-NEYMAN INTERVAL
##
## When TissueFraction_GMWM is OUTSIDE the interval [0.58, 1.13], the slope of
## mri_age_y is p < .05.
##
## Note: The range of observed values of TissueFraction_GMWM is [0.32, 1.00]
## SIMPLE SLOPES ANALYSIS
##
## Slope of mri_age_y when TissueFraction_GMWM = 0.5344591 (- 1 SD):
##
## Est. S.E. t val. p
## ------ ------ -------- ------
## 0.10 0.04 2.59 0.01
##
## Slope of mri_age_y when TissueFraction_GMWM = 0.6423103 (Mean):
##
## Est. S.E. t val. p
## ------ ------ -------- ------
## 0.03 0.04 0.83 0.41
##
## Slope of mri_age_y when TissueFraction_GMWM = 0.7501615 (+ 1 SD):
##
## Est. S.E. t val. p
## ------- ------ -------- ------
## -0.04 0.06 -0.74 0.46#save Johnson-Neyman plot
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/lag_ins_age_tf_int_JNplot.png")## Saving 7 x 5 in image#create interaction plot (make it pretty later)
#basic plot
Ins_age_tf_int_plot <- interact_plot(lag_ins_age_TisFrac_int, pred = mri_age_y, modx = TissueFraction_GMWM)
Ins_age_tf_int_plot
#linearity check
Ins_age_tf_int_plot_lin <- interact_plot(lag_ins_age_TisFrac_int, pred = mri_age_y, modx = TissueFraction_GMWM, modx.values="terciles", plot.points=TRUE, linearity.check=TRUE)## Medians of each tercile of TissueFraction_GMWM are 0.539, 0.655, 0.747Ins_age_tf_int_plot_lin## `geom_smooth()` using formula 'y ~ x'
##Print model summaries Print formatted model summaries and tables using report and sjplot package
library(report)
#library(sjplot) #try for pretty summary tables
summary(lag_naa_age_TisFrac)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 858.0 876.8 -424.0 848.0 313
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9925 -0.5718 -0.0199 0.5594 3.9440
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.1630 0.4037
## Residual 0.7184 0.8476
## Number of obs: 318, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 11.22872 0.43459 317.97211 25.838 < 2e-16 ***
## mri_age_y 0.17233 0.02961 317.97750 5.819 1.44e-08 ***
## TissueFraction_GMWM 3.35162 0.51255 313.38178 6.539 2.51e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.703
## TssFrc_GMWM -0.915 0.394report(lag_naa_age_TisFrac)## 'interpret_d()' is now deprecated. Please use 'interpret_cohens_d()'.
## 'interpret_d()' is now deprecated. Please use 'interpret_cohens_d()'.## We fitted a linear mixed model (estimated using ML and nloptwrap optimizer) to predict NAA_conc_molal with mri_age_y and TissueFraction_GMWM (formula: NAA_conc_molal ~ mri_age_y + TissueFraction_GMWM). The model included subj_id as random effect (formula: ~1 | subj_id). The model's total explanatory power is substantial (conditional R2 = 0.30) and the part related to the fixed effects alone (marginal R2) is of 0.15. The model's intercept, corresponding to mri_age_y = 0 and TissueFraction_GMWM = 0, is at 11.23 (95% CI [10.37, 12.08], t(313) = 25.84, p < .001). Within this model:
##
## - The effect of mri age y is statistically significant and positive (beta = 0.17, 95% CI [0.11, 0.23], t(313) = 5.82, p < .001; Std. beta = 0.32, 95% CI [0.21, 0.43])
## - The effect of TissueFraction GMWM is statistically significant and positive (beta = 3.35, 95% CI [2.34, 4.36], t(313) = 6.54, p < .001; Std. beta = 0.35, 95% CI [0.25, 0.46])
##
## Standardized parameters were obtained by fitting the model on a standardized version of the dataset. 95% Confidence Intervals (CIs) and p-values were computed usingreport_table(lag_naa_age_TisFrac)## 'interpret_d()' is now deprecated. Please use 'interpret_cohens_d()'.## Parameter | Coefficient | 95% CI | t(313) | p | Effects | Group | Std. Coef. | Std. Coef. 95% CI | Fit
## -----------------------------------------------------------------------------------------------------------------------------------
## (Intercept) | 11.23 | [10.37, 12.08] | 25.84 | < .001 | fixed | | 0.02 | [-0.11, 0.14] |
## mri age y | 0.17 | [ 0.11, 0.23] | 5.82 | < .001 | fixed | | 0.32 | [ 0.21, 0.43] |
## TissueFraction GMWM | 3.35 | [ 2.34, 4.36] | 6.54 | < .001 | fixed | | 0.35 | [ 0.25, 0.46] |
## | 0.40 | | | | random | subj_id | | |
## | 0.85 | | | | random | Residual | | |
## | | | | | | | | |
## AIC | | | | | | | | | 857.96
## BIC | | | | | | | | | 876.77
## R2 (conditional) | | | | | | | | | 0.30
## R2 (marginal) | | | | | | | | | 0.15
## Sigma | | | | | | | | | 0.85summary(lag_cre_age_TisFrac)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM + (1 | subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 760.2 779.0 -375.1 750.2 315
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.7493 -0.5056 0.0270 0.4419 4.8839
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.08547 0.2923
## Residual 0.54120 0.7357
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 7.04677 0.36923 319.70413 19.085 <2e-16 ***
## mri_age_y 0.05661 0.02497 319.62390 2.267 0.0241 *
## TissueFraction_GMWM 4.47115 0.43745 317.80861 10.221 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y
## mri_age_y -0.701
## TssFrc_GMWM -0.918 0.393report(lag_cre_age_TisFrac)## 'interpret_d()' is now deprecated. Please use 'interpret_cohens_d()'.
## 'interpret_d()' is now deprecated. Please use 'interpret_cohens_d()'.## We fitted a linear mixed model (estimated using ML and nloptwrap optimizer) to predict Cre_PCr_conc_molal with mri_age_y and TissueFraction_GMWM (formula: Cre_PCr_conc_molal ~ mri_age_y + TissueFraction_GMWM). The model included subj_id as random effect (formula: ~1 | subj_id). The model's total explanatory power is substantial (conditional R2 = 0.35) and the part related to the fixed effects alone (marginal R2) is of 0.25. The model's intercept, corresponding to mri_age_y = 0 and TissueFraction_GMWM = 0, is at 7.05 (95% CI [6.32, 7.77], t(315) = 19.08, p < .001). Within this model:
##
## - The effect of mri age y is statistically significant and positive (beta = 0.06, 95% CI [7.48e-03, 0.11], t(315) = 2.27, p = 0.024; Std. beta = 0.12, 95% CI [0.02, 0.22])
## - The effect of TissueFraction GMWM is statistically significant and positive (beta = 4.47, 95% CI [3.61, 5.33], t(315) = 10.22, p < .001; Std. beta = 0.53, 95% CI [0.43, 0.63])
##
## Standardized parameters were obtained by fitting the model on a standardized version of the dataset. 95% Confidence Intervals (CIs) and p-values were computed usingsummary(lag_cho_age_TisFrac_int)## Linear mixed model fit by maximum likelihood . t-tests use Satterthwaite's
## method [lmerModLmerTest]
## Formula: Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM + (1 |
## subj_id)
## Data: lag_dat
##
## AIC BIC logLik deviance df.resid
## 23.3 46.0 -5.7 11.3 314
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5480 -0.6259 -0.0670 0.5067 4.4339
##
## Random effects:
## Groups Name Variance Std.Dev.
## subj_id (Intercept) 0.01168 0.1081
## Residual 0.05179 0.2276
## Number of obs: 320, groups: subj_id, 95
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.45851 0.26231 312.15748 9.373 <2e-16
## mri_age_y 0.03027 0.03614 303.65993 0.838 0.403
## TissueFraction_GMWM 0.01894 0.40778 309.31482 0.046 0.963
## mri_age_y:TissueFraction_GMWM -0.11820 0.05899 301.29296 -2.004 0.046
##
## (Intercept) ***
## mri_age_y
## TissueFraction_GMWM
## mri_age_y:TissueFraction_GMWM *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) mr_g_y TF_GMW
## mri_age_y -0.943
## TssFrc_GMWM -0.981 0.948
## m__:TF_GMWM 0.896 -0.976 -0.942report(lag_cho_age_TisFrac_int)## 'interpret_d()' is now deprecated. Please use 'interpret_cohens_d()'.
## 'interpret_d()' is now deprecated. Please use 'interpret_cohens_d()'.## We fitted a linear mixed model (estimated using ML and nloptwrap optimizer) to predict Cho_GPC_PCh_conc_molal with mri_age_y and TissueFraction_GMWM (formula: Cho_GPC_PCh_conc_molal ~ mri_age_y * TissueFraction_GMWM). The model included subj_id as random effect (formula: ~1 | subj_id). The model's total explanatory power is substantial (conditional R2 = 0.28) and the part related to the fixed effects alone (marginal R2) is of 0.12. The model's intercept, corresponding to mri_age_y = 0 and TissueFraction_GMWM = 0, is at 2.46 (95% CI [1.94, 2.97], t(314) = 9.37, p < .001). Within this model:
##
## - The effect of mri age y is statistically non-significant and positive (beta = 0.03, 95% CI [-0.04, 0.10], t(314) = 0.84, p = 0.403; Std. beta = -0.32, 95% CI [-0.44, -0.21])
## - The effect of TissueFraction GMWM is statistically non-significant and positive (beta = 0.02, 95% CI [-0.78, 0.82], t(314) = 0.05, p = 0.963; Std. beta = -0.27, 95% CI [-0.38, -0.16])
## - The interaction effect of TissueFraction GMWM on mri age y is statistically significant and negative (beta = -0.12, 95% CI [-0.23, -2.13e-03], t(314) = -2.00, p = 0.046; Std. beta = -0.09, 95% CI [-0.18, -1.62e-03])
##
## Standardized parameters were obtained by fitting the model on a standardized version of the dataset. 95% Confidence Intervals (CIs) and p-values were computed using##Plot results of models including tissue fraction Estimate metabolite residuals using remef Include individual slopes (estimated via linear models by subject) using geom_smooth() Include overall fit line using geom_abline
#LAG - NAA
#pull fixed effects for overall model
fixef.lag_naa_TisFrac<-fixef(lag_naa_age_TisFrac)
#subset dataframe with only subjects included in NAA analysis
lag_dat_plot<- lag_dat %>%
filter(! is.na(NAA_conc_molal))
#calculate partial effect of age on tNAA (remove effect of tissue fraction)
lag_naa_partial <- remef(lag_naa_age_TisFrac, fix = "TissueFraction_GMWM") #does NOT remove random effects variance
#create plot to show NAA ~ age while controlling for Tissue Fraction
plot_lag_naa_age_tf<-ggplot(lag_dat_plot, aes(x=mri_age_y, y=lag_naa_partial, na.exclude = TRUE, color=as.factor(female))) +
geom_point() +
geom_smooth(method = "lm", se = FALSE, aes(x=mri_age_y, y=lag_naa_partial, group=as.factor(subj_id)), size = .2, show.legend = FALSE) +
labs(title = "Left Temporo-Parietal", x = "Age (years)", y = "tNAA | Tissue Fraction") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.lag_naa_TisFrac[[1]], slope = fixef.lag_naa_TisFrac[[2]])) +
theme_classic()
plot_lag_naa_age_tf## `geom_smooth()` using formula 'y ~ x'
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/lag_naa_age_partialef.png")## Saving 7 x 5 in image
## `geom_smooth()` using formula 'y ~ x'#LAG - Cre
#pull fixed effects for overall model
fixef.lag_cre_age_TisFrac<-fixef(lag_cre_age_TisFrac)
#subset dataframe with only subjects included in NAA analysis
lag_dat_plot_Cre<- lag_dat %>%
filter(! is.na(Cre_PCr_conc_molal))
#calculate partial effect of age on Cre (remove effect of tissue fraction)
lag_cre_partial <- remef(lag_cre_age_TisFrac, fix = "TissueFraction_GMWM") #does NOT remove random effects variance
plot_lag_cre_age_tf<-ggplot(lag_dat_plot_Cre, aes(x=mri_age_y, y=lag_cre_partial, na.exclude = TRUE, color=as.factor(female))) +
geom_point() +
geom_smooth(method = "lm", se = FALSE, aes(x=mri_age_y, y=lag_cre_partial, group=as.factor(subj_id)), size = .2, show.legend = FALSE) +
labs(title = "Left Temporo-Parietal", x = "Age (years)", y = "tCre | Tissue Fraction") +
scale_color_viridis(name = "Sex", labels = c("Male", "Female"), discrete = TRUE, option = "H", begin = .2, end = .8) +
geom_abline(aes(intercept = fixef.lag_cre_age_TisFrac[[1]], slope = fixef.lag_cre_age_TisFrac[[2]])) +
theme_classic()
plot_lag_cre_age_tf## `geom_smooth()` using formula 'y ~ x'
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/lag_cre_age_partialef.png")## Saving 7 x 5 in image
## `geom_smooth()` using formula 'y ~ x'#LAG - Cho
Cho_age_tf_int_plot <- interact_plot(lag_cho_age_TisFrac_int, pred = mri_age_y, modx = TissueFraction_GMWM, data=lag_dat, plot.points=TRUE, legend.main="Tissue Fraction GM/(GM+WM)", modx.labels=c(".535 (-1 SD)", ".642 (Mean)", ".750 (+1 SD)")) +
labs(title = "Left Temporo-Parietal", x = "Age (years)", y = "tCho (molal)") +
scale_color_viridis(option="mako", name="Tissue Fraction GM/(GM+WM)", begin= .32, end=1) +
theme_classic()## Scale for 'colour' is already present. Adding another scale for 'colour',
## which will replace the existing scale.Cho_age_tf_int_plot## Warning: Removed 1 rows containing missing values (geom_point).
ggsave("/Users/meaghan/OneDrive - University of Calgary/Preschool_data/MRS/results/lag_cho_age_tf_int.png")## Saving 7 x 5 in image## Warning: Removed 1 rows containing missing values (geom_point).